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ABSTRACT 


Ocean acoustic tomography permits the mapping of various properties of a body 
of water through indirect means. The technique utilizes travel time variations for an 
acoustic signal to determine the structure of the ocean medium via inverse mathematical 
methods. The scale of any tomography experiment is fundamentally limited by the 
signal to noise ratio at the receiver. Through the use of a near vertical acoustic array, 
normal mode modeling of the local environment and a modal beamformer, array gains 
are possible which greatly extend the maximum separation between source and receiver. 
Additionally, the technique provides temporal resolution of the modal components of 
the arriving signal. 

A time domain modal beamformer for a near vertical acoustic array has been 
developed. It has realized a noitnal array gain of 6 dB for the Heard Island 
Experiment vertical array deployed off California. The primary obstacle to the 


technique remains inadequate array geometry description. 
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I. INTRODUCTION 


A. THESIS SUMMARY 

The objective of this thesis is to develop a software 
package to beamform acoustic signals used in ocean acoustic 
tomography. The goal of tomography signal processing is the 
precise measurement of acoustic travel time. Sound speed is 
a well understood function of temperature, pressure and 
Salinity. Therefore, a fluctuation in acoustic travel time is 
indicative of changes in the environment through which the 
acoustic energy has passed. Once the travel time fluctuations 
are known, inverse mathematical techniques can be used to 
Miner Various properties of the ocean medium. [Ref. 1) 

One source of travel time measurements is through the use 
of explosive devices. Although such signals are easy to 
generate, they may not provide the required precision owing to 
dispersion of the various frequency components in the 
impulsive signal. Additionally, such signals are difficult or 
HmpOosslplewto replicate. [Ref. 2] 

A better method that has been employed in recent years is 
the use of maximal-length sequences. The technique utilizes 
a pseudorandom, binary sequence as the phase modulation for an 
electronically generated bandpass signal. Maximal-length 


sequences are well suited to travel time measurements by 


virtue of their detereminiret ire merctucc, autocorrelation 
properties and simplicity. 

The goal of the work described in the following pages is 
to provide sufficient gain to permit the recovery and decoding 
necessary for accurate travel time measurements for acoustic 
paths of extreme length (in this case 10,000 nautical miles). 
Specifically, the programming package should be able to: 

@ Exploit the spatial structure of the incident wave field 
to discriminate among signals arriving from different 


directions. 


@e Utilize information regarding time varying array tilt and 
depth to estimate the array geometry. 


e Provide a stable virtual array by uSing the array geometry 
and modal structure of the immediate environment. 


e Distinguish among the various modal components of the 
target signal. 


@e Provide both time domain and frequency domain analysis. 

@e Remain sufficiently flexible so as to permit the 
processing of arrays of differing construction in future 
tomography experiments. 

@ Provide some measure of fault tolerance with respect to 
the receiving array. 


The remainder of this thesis is structured as follows: 


e Chapter II. Acoustic Wave Propagation Theory 
e Chapter III. Modal Beamforming 
e Chapter IV. The Heard Island Array 


e Chapter V. Results and Conclusions 


Chapter II presents an introduction to acoustic wave 
propagation. The approach taken is that of normal mode 
propagation in a range independent channel. 

The third chapter is an introduction to the concept of 
modal beamforming. It reviews frequency domain and time 
domain beamformers. The algorithm utilized in this study is 
derived. 

Chapter IV describes the construction and subsequent 
deployment of the receiving array used in this experiment. 
Data acquisition and preprocessing of the acoustic signal are 
discussed. 

The final chapter presents a description of the array 
dynamics encountered and frequency domain output from the 
software developed. Additionally, proposals for system 
integration and performance improvements are detailed. 

Appendix A contains the source code which has served as 
the body of this work. In addition to the beamformer, various 
utility programs are included which provide for array geometry 


data reduction. 


B. THE HEARD ISLAND EXPERIMENT 

The performance of the software package is evaluated on a 
raw data set (vice synthetic data). The data was acquired 
during the Heard Island Experiment which occurred during the 
period January 23 - February 2, 1991. The purpose of the 


experiment was to determine the reliability of global acoustic 


paths for tomographic analysis. Specifically, it is desired 
to establish the viability of these transmission paths for a 
proposed multi-national attempt to detect a decreasing trend 
in acoustic travel time. such a change would indicate an 
overall warming along the path, providing the first convincing 
evidence to the existence of glebal warnings) Ret. iE 

The signals were transmitted from the vicinity of Heard 
Island in the southern Indian Ocean. The site was selected 
because it is central to a web of open water paths extending 


through all the worlds oceans as shown in Figure 1.1. 
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Figure 1.1: Heard Island Raypaths 


The transmitting ship was the R/V Cory Chouest. The 
transmitter utilized a ten element vertical array with a 
maximum of five elements active at any time. A nominal source 
level of 209 dB re 1 upPa at 1 meter was realized. Individual 


projectors employed were HLF4LL very low frequency sources. 


A variety of signal types were sent, including continuous 
wave and maximal-length sequences. A carrier frequency of 57 
Hz was chosen both for its low absorption losses and the 
ability to distinguish it from the 50 and 60 Hz frequencies 
generated by power plants world wide. 

The Monterey component of the Heard Island Experiment 
involved a collaboration between the Naval Postgraduate 
School, the Monterey Bay Aquarium Research Institute, the 
Massachusetts Institute Gt Technology, Woods Hole 
Oceanographic Institution, SAIC and the Moss Landing Marine 
Laboratory. The R/V Point Sur deployed a 32 element vertical 
array approximately 70 nautical miles south-west of Monterey 


Bay. 


II. ACOUSTIC WAVE PROPAGATION THEORY 


A. ACOUSTIC PROPAGATION IN AN OCEAN WAVEGUIDE 
1. The Inhomogeneous Wave Equation 


The homogeneous linear wave equation [Ref. 3], 





nO 
a ee ' (2.1) 
Co Oe 
incorporates two implicit assumptions. The first restricts 


the sound speed to a constant value in the region of interest. 
Following the development of Coppens [Ref. 4], if the sound 
speed is now permitted to vary in space (but not time), then 
c=c(x,y,z) may be substituted into the wave equation without 
loss of generality. 

The second implicit assumption restricts the use of 
the wave equation to sound fields which are free of sources. 
Before including a source termpeene Selection of a Speciree 
coordinate system is appropriate. 

Given that sound in the ocean does not’ spread 
spherically ina free field, but radially with upper and lower 
boundaries, implies the use of a cylindrical coordinate 
system. To further simplify the problem to one sufficient for 
this work, the sound field shall be assumed to exhibit radial 
symmetry about the source, and sound speed variations shall be 


TEStricted to (Gepen mmciuc Based on these restrictions, 


application of the Lapacian yields the homogeneous, range 


independent wave equation, 





iO} ..0 Se 1 & 
=a a)? ae? c2(z) ate 


The inclusion of a source term requires that the 
inhomogeneous wave equation reduce to homogeneous form in 
regions which are free of sources of acoustic energy. This 
requirement implies that the delta function for a point source 


Beeated at r-0 and z=z, have the form [Ref. 3] 


S(r-r,) = —~—8(r) 8(z-z,), (2.3) 


PAs ie 





where r is the spherical radius vector, r is the cylindrical 
radial coordinate and z is the cylindrical depth coordinate. 

Inclusion of this term in the wave equation yields the 
Helmholtz equation for a monofrequency point source of unity 


amplitude, 
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(2.4) 

2. Solutions to the Inhomogeneous Wave Equation 
Having adopted a cylindrical coordinate system, 
equation (2.4) can be solved by separation of variables. 


Furthermore, the complete solution can be treated as a linear 


combination of normal modes, 
p(r, 2, t) se)? Gee) (255 
Mm 


The normal modes 9(2Z_) 9fonm—@panemerelnonormal | setCumes 
eigenfunctions in z which satisfy the sourceless Helmholtz 


equation, 








2 
q ae w? he Zm = O, (2.6) 
az Caz 


subject to the appropriate surface and bottom boundary 


conditions, and the orthonormal condition, 


[2,.(2) 2,2) dz = 8 (2.7) 


nm! 


where k, is the eigenvalue for the m" eigenfunction (or normal 
mode) and oO. is™tne Kronecker Je) ea aumc eter, 

A closed form expression for the R_(z) may be obtammea 
by substitution of equation (2.6) into the inhomogeneous wave 
equation (2.4), multiplication of all terms by 7 eimvegraeeag 
over z and application of the normalization condition (2.7). 


The result, after manipulation, is 





ae ane 2 Z 2.8 
|| ie +KOR = -—§6(r)Z ; (2.8) 
ie =| on a i PEW) 
the inhomogeneous Bessel's equation. It has the known 
solution 
Ror) = =9 0 oe see (2.9) 


where H,“ (kr) is the Hankel function of the second kind and 
ema@er Zero. Substitution of this form into equation (2.5) 


yields 
Pitgeeee rome ) Zz) 2.) He) (Kr) . (2.10) 
Mm 


micccOmeeo Tm oOLULTONS fOr “the @2Z8(z) anes either 
difficult or impossible for all but select’ cases. 
Fortunately, efficient numerical algorithms exist which 
provide for the rapid solution of the depth dependent 


mumetions. [Ref. 5] 


III. MODAL BEAMFORMING 


A. PRELIMINARY CONCEPTS 

Propagating waves can be modeled as functions of both 
space and time. Consequently, the attributes of such models 
can be used to extract information from real wave fields. 
Beamforming exploits the temporal and spatial characteristics 
of a specified environment to enhance a particular aspect of 
the wave field while attenuating undesirable components. Such 
an operation can be loosely termed constructive reenforcement 
of the desired signal or noise rejection. The information 
exploited in this development includes the modal structure of 
the immediate environment and the mnstrument data supplicase 


the array itself. 


B. TIME DOMAIN BEAMFORMING 
The acoustic signal received at a given hydrophone can be 


expressed as 


D(C, 2) = "SCS ea ey ee (3.1) 


where p is the pressure at the hydrophone, S is the desired 
Signal and N is the local noise field. All are functions of 
time and hydrophone location. 

If one now considers the signal received at an array of 


hydrophones, equation (3.1) takes the form 
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ie oC) Ne); (3.2) 


where p, S and N are all N X 1 vectors representative of the 
individual array elements. The location of individual 
elements is implied by the vector index for a given 
hydrophone. 

Classical (or plane wave) beamforming incorporates the use 
of a complex steering vector (E) to impose a phase shift on 
feneeindividual elements in order to enhance the sensitivity of 


the array to signals propagating in a specific direction, 


E=| |, (3.3) 








The output signal from such an array is 


ie (ate Me eae Sua (te) (3.4) 


The subscript (pw) indicates that the beamformer assumes the 
presence of plane waves. Additionally, the implicit 
assumption is made in the above expression that the 
autocorrelation function of the noise field will be 
identically zero for any two (distinct) hydrophones. The 
Success or failure of the beamformer will be a direct 


reflection of the validity of these assumptions. 


ie 


If the desired signal is not monofrequency, but of finite 
bandwidth, simple phase shifting may not be suitable. This is 
particularly true if the target signal is carrying information 
such that the communication bit time is comparable to, or 
greater than the time taken for the wave to propagate across 
the array. Under such circumstances, the beamformer would 
impose phase shifts in distinct communication bit segments 
Simultaneously. The resulting distortion to the bandpass 
Signal is unacceptable in most communications applications. 
Beamforming under these conditions requires the application of 
true time delays vice phase shifts. 

Assuming an array of spatial extent such that time domain 


beamforming is required implies a summation of the form 


GS Gz 0) sige) le, (3.5) 


M= 


ae i 
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n= 


where n is the index used to describe the array elements and 


h 


t. represents the steering delay applied to the n‘ element 


rn 
(Ref. 6). This arrangement places no additional restrictions 
on the target signal characteristics or array constructions 


a given physical system. 


c: THE MODAL BEAMFORMER 
The signal model employed above is incomplete, in that it 
does not account for the nonuniform distribution of acoustic 


energy in the sound channel. Vertical arrays inherently 
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sample the modal structure of the immediate environment. 
Therefore, one needs to take the information provided by 
normal mode theory into account when the acoustic array spans 
a region of nonuniform energy density. 

The appropriate refinement to the complex steering vector 


which incorporates the modal nature of the signal arrival is 


ARGZ,) e! Fu Ze a) e! Viz , : : Zz) el Yim 
Le ca eo V2: Ze (z,) et Vz ae Sat on Ze Zs) el Yam 
oes | , eee , (3.6) 
. 1 
eee ee en (Zp) er tele. we ZZ) Ce 
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where the mt? column of the matrix is the generalized steering 


a 


wector for the m, normal mode. The individual steering 


h 
feelous Now inevude amplitude factors (2 (z.)) that reflect the 
value of the eigenfunction at the various hydrophone depths 
(Ref. 7]. The nt row contains the values of those steering 
vectors at the depth of the nt hydrophone. The phase terms 
meeeccount £or phase Shifts required for the beam steering 
applicable to all modes (as required by array geometry) and 
phase delays among the various modes due to their individual 


propagation speeds. Specifically, the ¥, take the form 


Won = Ger a ; (3.7) 
where 8, is determined by the hydrophone location anda, is 


determined by the differences in propagation speeds over the 


transmission path for the normal modes present. 


3 


Modifying the steering matrix to implement time delays 
(vice phase shifts) is straight forward. A beamformer could 
be implemented by summing the output of individual hydrophones 
(p(t)) at this point. However, to do so implies a 
deterministic evaluation of thevame Such anvevaluation woud 
require near perfect knowledge of the conditions over the 
transmission path prior to commencement of the experiment. 
Given that the aim of tomography is the. mapping of these 
conditions, such an approach is inappropriate. 

An alternate approach is to utilize the columns of the 
matrix individually to discriminate among the modal components 
of the signal. Since all components of a given modal steering 
vector share a common value of a,, these phase terms may be 
discarded from the individual elements of the steering vector. 
As a result of this simplification, the modal steering vector 
is completely determined by the Z (obtained from normal mode 
modeling of the immediate environment) and the desired look 
direction (determined by array geometry). Having established 
the methodology, and replacing the phase weights with steering 


delays, the output of the modal beamformer is 
N 
DE Gey =e 2 Cz, noma ele (3.8) 
net 


where the subscript (m) indicates that the desired output is 


h 


the component of signal representing the m‘” normal mode. 


Expanding this form to explicitly “inclide thew tem- 
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representing the signal and noise components of the acoustic 


field yields 


N N 
Die eZ 2) Sa(G=%.) +) 2.( 2.) Ne tae. (3.9) 
n= 


n#1 


where S, and N, represent the amplitudes of signal and noise 
at the n™ hydrophone. The delay and weighting of the noise 
received at the hydrophones has no undesirable effects. If 
the noise between any two hydrophones was uncorrelated prior 
to the operation, then it will remain so after weighting. 
Again, the degree to which the ambient noise field is 
uncorrelated will ultimately reflect upon the performance of 


the beamformer. 


in 


IV. THE HEARD ISLAND VERTICAL ARRAY 


A. ARRAY CONSTRUCTION 

The array deployed for this experiment consisted (ome 
hydrophones deployed vertically, each having a nominal 
sensitivity of -170 dB re 1 V/ypPa. The element spacing was 45 
meters, with hydrophone number one occupying a design depth of 
345 meters. The nominal design depth for hydrophone number 32 
was 1740 meters. 

Instrumentation on the array consisted of two sensors. 
The upper sensor was located 4.0 meters above hydrophone 
number one. It recorded ambient pressure, temperature, tilt 
and current velocity. The lower sensor was located 5.0 meters 
below hydrophone number 20. It recorded tilt (ancludime 
direction), pressure, temperature and conductivity. This 
sensor appears to have suffered a casualty during array 
deployment which rendered its tilt data suspect. 

Flotation was provided by one main syntactic float and 28 
Syntactic football floats. The design called for the main 
float to reside at a nominal depth of 230 meters. 
Approximately half of the football floats were submerged. 
This arrangement rendered the array neutrally buoyant, thus 
providing a meaSure of isolation from surface wave effects. 


The array was directly tethered to the R/V Point Sur with 
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floatation devices on the surface portion. Each float was 
equipped with a flasher, a flag and a radar reflector. The 
float nearest the ship had the addition of a radio beacon. 
The procedure called for the ship to remain approximately 1200 
meters from the position immediately above the array. During 
the experiment the R/V Point Sur was configured to remain 
quiet and dead in the water, except as required to keep clear 
of the array. Figure (3.1) illustrates the general 


configuration. 
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Figure 3.1: Heard Island Receiving Array 


B. INSTRUMENT DATA 

Ideally the receiving array, when deployed, is oriented 
vertically. However, due to currents, the requirement to 
maintain control of the research vessel and other factors 
beyond the control of the research team, the array is often 
tilted from the vertical. Given the length of the ammay 7a 
1.0° tilt results in a 24 meter horizontal displacement 
between hydrophones one and 32. This displacement is nearly 
one (carrier frequency) wavelength. Additionally, the 
direction in which the array tilts affects the component of 
horizontal displacement which is collinear to the propagation 
direction. Horizontal displacements perpendicular to the 
direction of propagation are acceptable if one assumes that 
the wave field is fairly wunitorm 1s thicw@ic@eecerere 
Therefore, the array steering utilized in this work seeks only 
to correct for the component of displacement in the direction 
of signal propagation. 

The speed with which the wave field propagates is required 
prior to the calculation of steering delays. The appropriate 
speed to use in this case is the group speed of the carrier 


frequency for the mode in question, defined as 


Com = Se (3.10) 


m 


h 


where Co, 15 the local group speed for the m‘” normal mode, o 


is the radian frequency and x, is the appropriate eigenvalue. 
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The amplitude weight applied to a given hydrophone is a 
function of the mode number sought and hydrophone depth. The 
Meade shapes (25(2)) are required prior to beamforming. The 
numerical algorithm used to obtain the eigenvalues and 
Seegentunctions for this work was developed by Chiu and Ehret 
(Ref. 5]. It utilizes finite differencing of the propagation 
vectors for a given frequency and sound speed profile to 
calculate the modal structure of the local environment. The 
sound speed profiles were measured during the conduct of the 
experiment by the research team aboard the R/V Point Sur. 
Figure 3.2 is the sound speed profile used for this analysis. 
Figures 3.3 through 3.5 are the output of the normal mode 
model. Amplitude weights for the array elements were 


calculated from these data sets. 


Ju 


sound Speed Profile 
Heard Island Receiving Array 





Depth (meters) 


-2000 
1475 1480 1485 1490 1485 1500 


Sound Speed (meters/second) 


Figure 3.2: Sound speed profile at receiving array. 
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Figure 3.3: Normal modes one and two. 
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Figure 3.4: Normal modes three and four. 
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Mode Five Mode Six 
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Figure 3.5: Normal modes five and six. 


C. DATA ACQUISITION AND PREPROCESSING 

Once deployed, the receiving array recorded numerous data 
sets, each lasting in excess of 60 minutes. Each recording 
comprises 32 channels of acoustic data. The hydrophone output 
was ported to a patch panel and preamplifier where a fixed 
gain of 26 dB was applied. Following this, the signal passed 
through a variable gain amplifier which automatically set the 
Signal gain to optimize the 10 volt dynamic range of the 
analog to digital converter. This process had the duel ebtece 
of enhancing the precision of the data and preventing 
Saturation of the A/D converters. Subsequent normalization of 
the data set to hydrophone output levels was accomplished by 
an implementation of a program written by Keith Von der Heydt 
of Woods Hole Oceanographic Institution. 

Analog to digital conversion utilized a sampling frequency 
of 228 Hz. The signal was passed through an analog bandpass 
filter prior to conversion. This prevented aliasing of 
frequencies above Nyquist into the data set. Low frequency 
noise components were also attenuated. The cutoff frequencies 
for the bandpass filter were 10 and 80 Hz. Data storage 
utilized both optical disk and magnetic tape. Data 
acquisition and preprocessing of the data sets used in this 
work was conducted by Miller, Chiu, Frogner and others. 


Figure 3.6 illustrates the equipment alignment. [Ref. 8] 
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Figure 3.6: Data acquisition system 
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V. RESULTS AND CONCLUSIONS 


The desired result of this research is the development of 
a software package capable of beamforming a tomography signal. 
The beamformer should discriminate among the various modal 
components of the signal, thus enabling one to treat the modes 
individually. Due to the inherent nature of vertical arrays, 
the software needs to incorporate information pertaining to 
the time varying array tilt and individual hydrophone depths. 
To this end, the beamformer utilizes a duty cycle. Each en 
cycle commences with an assessment of array geometry. This 
information is used to determine the appropriate steering 
delays and hydrophone weighting coefficients. Once 
established, the beamformer processes 60 seconds of acoustic 
data before repeating the procedure. In this manner, sen 
software package implements a gquasi-stable virtual array in 
which the virtual aperture remains fixed in space. 

The results presented examine two performance aspects. 
The first pertains to array geometry description. Calculation 
of the steering delays and amplitude weights is detailed. 
Additionally, the acoustic output of the array is compared to 
that of a single hydrophone on the sound channel axis, thus 
allowing for a quantitative evaluation of the ultimate 


performance of the beamformer. 
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A. HYDROPHONE DELAY AND WEIGHTING 

Information regarding array geometry was provided by two 
instrument packages. The upper sensor was located 4.0 meters 
above hydrophone number one, the lower sensor 5.0 meters below 
hydrophone number 20. Figure 5.1 illustrates the instrument 
configuration and coordinate system used in determining the 


array's spatial orientation. 
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Figure 5.1: Array Instrumentation 
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1. The Upper Sensor 

Information provided by the upper instrument includes 
the component of current velocity in both the north-south and 
east-west directions, ambient pressure, temperature and tilt. 
The direction in which the array tilts is determined by the 
current direction. The primary assumptions being that the 
current acting on the array is independent of depth and that 
the direction in which the array tilts is completely 
determined by the current direction. 

In making these assumptions, one considers the primary 
source of current measured by the instrument to be a result of 
ship positioning adjustments. As mentioned previously, the 
R/V Point Sur was required to periodically move the ship in 
order to remain clear of the array. Such movements were 
intended to straighten the tether between the array and the 
towing rig, thus preventing the tether from becoming fouled in 
the ship's propellers. These adjustments resulted in a 
tensioning of the tether. It is this tensitonings that pre 
considered to have resulted in array movement and consequent 
water flow across the instrument package. Based on this 
assessment, the assumption that tilt direction was determined 
by the flow measured at the upper instrument is consistent 
with the conditions of the experiment. Higher order models 
could be derived for array dynamics, however this research has 


utilized the first order model of a linear array. 
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The next variable pertinent to array geometry is the 
pressure measured by the upper sensor. From this it is 
straight forward to calculate the sensor's depth using the 


following relationship, 


P-latmosphere (5.1) 
Pg 


Z2nt 


where P is the pressure “(in pascals) measured by the 
instrument, p is the density of seawater and g is the force of 
gravity. 

Figures 5.2 and 5.3 present the data received from the 
Meper instrument package during two of the acoustic recording 
periods. iWtemm@cecoraing beginning 00:05 (GMT) January 27 
coincides with the reception of a m-Sequence of length 2047. 
Pepeeneinuous wave signal was received during the recording 
beginning 15:05 (GMT) January 27. They are representative of 
the most and least favorable conditions encountered during the 


course of the experiment. 
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igure 5.2: Array geometry for 00:05 (GMT) January 27. 
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2. The Lower Sensor 

The lower sensor was intended to provide a redundant 
source of array positioning data. Specifically, it measured 
tilt in the north-south and east-west directions, temperature, 
pressure and conductivity. As discussed previously, the lower 
instrument appears to have suffered a casualty which rendered 
the tilt data suspect. 

Figure 5.4 indicates that an event occurred which 
caused the sensor to provide data which is not consistent with 
the assumed array geometry. One scenario explaining the event 
is that the sensor failed shortly after entering the water. 
Another possible explanation is that the array became 
entangled in itself upon descent. 

All array tilt data was stored for subsequent analysis 
ashore. Therefore, the research team had no knowledge of the 
apparent failure during the conduct of the experiment. As 
there was no reason to suspect the instrument, no 
investigation was made to determine the source of the 


inconsistencies. 
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Figure 5.4: 


Lower sensor output upon deployment. 
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Depth (meters) 


3. Steering Delays 

Time delays must be applied to each of the hydrophones 
in order to correct ion  auray see However, only the 
component of horizontal displacement collinear to the signal 
propagation direction is used to calculate the time delay. 
This requirement mandates apriori knowledge of the propagation 
direction. Prior work has determined that the expected 
direction from which the signals would arrive is 217.0° true 
(Ref. 9]. 

In addition to array geometry and propagation 
direction, one requires a nominal propagation speed to convert 
the horizontal displacements to time delays. The appropriate 
speed to use is the mode group speed, as previously discussed. 
Group speeds for this study were obtained from eigenvalues 
calculated by the normal mode model courtesy of Chiu and Ehret 
(Ref. 5]. The mode group speed of the carrier frequency was 
used to calculate time delays. Figure 5.5 presents the mode 
group speeds as a function of frequency for the lowest five 


modes. 


34 


Mode Group Speeds 
1482.0 


1481.5 : ae 
1481.0 : i 
1480.5 


1480.0 


Group Speed (meters/sec) 


1479.5 


Mode Five 





1479.0 
0 0 0 DM 0 0 OH 


Frequency (Hz} 


Figure 5.5: Mode group speeds. 


4. Hydrophone Amplitude Weighting 
Application of hydrophone weights utilizes the array 
geometry to determine the depth of each array element. The 
value of the depth dependent function’(Z_) is then assigned to 
the appropriate hydrophone output. This manner of weighting 


" normal mode. 


enhances the array sensitivity to the m* 

An additional consideration encountered in hydrophone 
weighting is fault tolerance. Vertical acoustic arrays 
operate ina hostile environment, thus individual hydrophones 
may fail. The amplitude weights applied to the array must 
cope with the failure of individual elements. To this end, 
the beamformer applies a weight of zero to elements which have 
been evaluated as unreliable. 

Log books kept by the research team aboard the R/V 
Point Sur indicating suspect hydrophones were verified by post 
processing statistical analysis of the individual elements. 
Several hydrophones were found to have failed during the 
course of the experiment. Specifically elements 18 through 28 
were evaluated as _ unreliable. As a result of this 
determination, only hydrophones one through 17 were used to 


beamform the target signals. 


B. ACOUSTIC PERFORMANCE 
The preliminary data analysis seeks to quantify the 
beamformer's performance with respect to a single hydrophone 


located on the sound channel axis. The power spectrum 
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Calculated from the output of hydrophone five was compared to 
the spectra calculated from the beamformed output of channels 
one through 17. The array gain is defined as 


SNR 
AG = 10109,0| Zr" (5.2) 
Pe 


where SNR, is the signal to noise ratio of the beamformed 
output and SNR, is the signal to noise ratio of the n‘™ 
hydrophone [Ref. 10]. The SNR for each case was calculated by 
integration of the power spectra over a 6.0 Hz bandwidth 
centered on the continuous wave frequency. 

The array gain was calculated for numerous power spectra, 
each representing 60 seconds of acoustic data. Additionally, 
the gain was calculated for the lowest five normal modes 
present. Values for array gain range from a low of 2.5 dB to 
a high of 8.5 dB. The nominal array gain is 6 QB. 

Reasons for the range in values include temporal 
variability in the modal components of the signal, errors 
introduced in the steering delays and errors in amplitude 
weighting. Time delay and amplitude weight errors are the 
direct result of assumptions made regarding array geometry. 
The primary source of geometric error is considered to involve 
the assumption that tilt direction is completely determined by 
current direction. Array tilt bearing may not be aligned well 


with current direction during periods when the direction of 


flow across the sensor is varying. Therefore, during periods 


oy 


of variable current flow, marginal performance can be expected 
of the beamformer. 

Figures 5.6 through 5.11 present representative power 
spectra of single channel data and beamformed data for modes 
one through five over the entire range of frequencies 


analyzed. 
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Figure 5.6: Hydrophone five power spectrum. 
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Figure 5.7: Mode one power spectrum. 
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Figure 5.8: Mode two power spectrum. 
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Mode Three Power Spectrum 
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Figure 5.9: Mode three power spectrum. 


Mode Four Power Spectrum 
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Figure 5.10: Mode four power spectrum. 
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Mode Five Power Spectrum 





5.0x10~ 
4.0x10~ ceceonntinetdiesseneeeieneee Sere becceeeeeneeeedecenaeeecnenee Daher 
ree Me I 


2.0x10% 


Power Spectral Density (V2 / Hz) 


HCI Gl Mens Seetacet cree arsacdees Boeaee-- beccondfeceeeesfueeeescecnneedaceeseeeeseee 


0.0x10° - ee piel Del ay RM a raha 
O 2O 40 60 So 1006 120 


Frequency (Hz) 


Figure 5.11: Mode five power spectrum. 


As shown in the previous figures, the noise field at the 
beamformer 1S not incoherent. This fact, coupled with 
imperfect geometric array description, detracts from the 
calculated array gain. 

One feature worth noting is the decreasing strength of the 
undesirable tonal components with increasing mode number, 
specifically those in the 40 to 50 Hz range. Presumably, 
these components were radiated by the R/V Point Sur. One 
hypothesis is that the strength of these signals in mode one 
is the result of grating lobes in the array's beam pattern for 
this mode. The amplitude weights applied to the hydrophones 
are similar to those for a plane wave beamformer, in that 


there is no phase shift applied along the array in the form of 
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coefficients of differing sign. Therefore, spatial aliasing 
is possible for frequencies with wavelengths less than twice 
the hydrophone spacing. For the Heard Island array, this 
includes frequencies ranging from 17 Hz to the bandpass cutoff 
Paequency of GO Hz 

As the mode number sought increases, these frequency 
components diminish in strength. This indicates that the 
rejection of coherent noise of local origin improves as the 
number of 180° phase changes applied to the hydrophone output 
increases. 

Also shown in these figures is the presence of distant 
sources of coherent noise. These other tonal components are 
seen to vary both as a function of mode number and time. This 
observed behavior is consistent with modal propagation in a 
range dependent ,environment. Merchant traffic in shallow 
water and offshore drilling activity are two possible 
candidates for sources of coherent noise which becomes coupled 
into the sound channel. 

Figures 5.12 through 5.17 present the same power spectra 
in the immediate vicinity of the carrier frequency (57 Hz). 
Again, the beamformed output indicates successful attenuation 
of undesirable components of a coherent noise field, 


presumably Of ,locaieonr1 gr. 
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Figure 5.12: Hydrophone five power spectrum. 
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Figure 5.13: Mode one power spectrum. 
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n Mode Two Power Spectrum 
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Figure 5.14: Mode two power spectrum. 
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Figure 5.15: Mode three power spectrum. 
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Mode Four Power Spectrum 
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Figure 5.16: Mode four power spectrum. 
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Figure 5.17: Mode five power spectrum. 
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C. CONCLUSIONS 
Modal beamforming is a valuable technique for use in long 
range tomography experiments. The signals presented in this 
work traversed the longest transmission path for which 
reception was attempted during the course of the Heard Island 
Experiment. Although no travel time estimates have been made, 
the successful detection of the signals indicates that large 
scale tomography experiments are feasible. 
The primary concern of this work was the development of 
a system capable of detecting the acoustic signals emanating 
from the vicinity of Heard Island. This has been 
accomplished. Additionally, the data describing the complex 
dynamics of the near vertical acoustic array has implications 
1or the design and construction of vertical arrays for use in 
future tomography experiments. The most troublesome aspect of 
this work has been the accurate determination of array 
geometry, specifically tilt direction. Describing tilt 
direction in terms of measured current at one point on the 


array is considered inadequate. 


D. RECOMMENDATIONS 

Improvements to this implementation include physical 
aspects of the array itself and performance of the numerical 
algorithm. Large vertical arrays are particularly sensitive 


to errors in estimated geometry. In this experiment, an error 
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eaeo.25° in tilt (assuming perfect knowledge of tilt 
direction) creates a phase error of 45 degrees between the 
most distant elements. This error is more destructive to the 
higher modes owing to their greater spatial extent. The 
problem is aggravated by questionable tilt direction data. 
Construction of future vertical arrays should incorporate 
multiple instrument packages which directly measure tilt and 
direction. As this was the case with the lower instrument on 
the Heard Island array, it is most unfortunate that it failed. 
Peeweionally, an investigation into the utility of higher 
order models describing array shape would be illuminating. 
Incorporating a quadratic model for array shape would not 
impose a significant computational load on the beamformer and 
may improve the estimated position of individual hydrophones. 
The current implementation of this beamformer utilizes 
third order polynomial interpolation during the application of 
time delays to individual hydrophone outputs. As such it 
represents the greatest computational load in the program. 
The interpolator is based on Neville's algorithm [Ref. 11}. 
No attempt was made to minimize execution time. As a result, 
the software requires approximately four hours to process one 
hour Of acoustic data. A more efficient interpolator may 


reduce execution time sufficiently to permit this technique to 
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be used for real time modal beamforming. To this end, the 
software is assembled from modular components, thus allowing 


for improvements on a function by function basis. 
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APPENDIX A 


The following programs were developed during the course of 
this research. They include the modal beamformer and 
associated utility programs. These programs may be obtained 
by contacting James H. Miller at his address in the initial 


Gistribution statement. 


A. THE MODAL BEAMFORMER 
1. Operational Considerations 
The modal beamformer was designed with portability a 
Peame consideration. As a result, there are numerous program 
parameters which must be set by the end user. Among these are 
the characteristics of the physical array and data acquisition 
system. All such user selectable parameters are contained 
solely within the program definitions. The parameters appear 
below in the same order in which they appear in the program. 
@ UNIX VERSION: selectively compiles either unix or ansi 
compatible code blocks. 
Wee eee Se leCeS Output an e€ither ascii or binary format. 
Ascii is useful for data export, binary saves storage 


space. 


@ SIGNAL: permits user to enable or disable the beamformer's 
time series output. 


@ SPECTRUM: permits the user to enable or disable the 
beamformer's power spectrum output. 
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LOWER SENSOR: directs the program to load or ignore lower 
tilt sensor data. This option should remain OFF as there 
is no reliable lower sensor data. 


VALIDATE: enables program validation. If selected the 
program will create output files which list the array 
geometry, steering delays, hydrophone weights’ and 
calculated mode group speed. This feature is useful for 
program debugging. 


ERROR_ESTIMATE: if selected, this feature causes the 
polynomial interpolator to store the upper bound on the 
mean squared error in the steering delays. This feature 
doubles the output Of tne beantormer. 


ON: logical switch for program control. 
OFF: logical switch for puogran corenel. 
INTERPOLATE: selects interpolator. This feature permits 
the installation of an improved interpolator without the 
requirement to search for and replace each call to the 


Funct vor. 


ORDER: indicates the degree of the polynomial used in the 
polynomial interpolator. 


STEP: indicates the number of points on either side of a 
sequence to be taken when estimating derivatives. 


TINY: prevents division by zero in floating point 
operations. 


PI: used for trigonometric recursions. 
RADIAN: used for degree to radian conversions. 


OFFSET: indicates the distance (in meters) between the 
number one hydrophone and the upper tilt sensor. 


DELTA _R: indicates array element spacing (in meters). 
CTD OFFSET: indicates the difference in depth between the 
sound speed profile's first data point and the depth 
increment used in the profile (in meters). 

SSP_LENGTH: indicates the maximum number of data points 


that an eigenfunction may contains Ine only Beeserictier 
applicable to this parameter is available memory. 
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@ EIGVAL_ LENGTH: indicates the maximum number of 


eigenvalues that the program will accept. The only 
restriction applicable to this parameter is available 
memory. 


@e LOOK DIRECTION: the compass direction (in degrees true) 
from which the signal arrives. 


Seer enUrrE Rs sthe burter length allocated to tilt data. 
The only restriction applicable to this parameter is 
available memory. 


e BUFFER TIME: the length of time (in seconds) represented 
by the acoustic data input buffer. This period represents 
the output buffer length plus one second on either end to 
permit steering delays and/or advances. The array is 
currently capable of end firing a 1500 meter acoustic 
array. This parameter must be an integer. Additionally, 
60 must be an integer multiple of (BUFFER _TIME-2). The 
array has a duty cycle of one minute between design 
changes. There must be an integer number of output 
buffers per duty cycle. 


@e F SAMPLE: sampling rate of the data acquisition system. 
This value must be an integer. 


@e CHANNELS: the number of hydrophones on the array. 

e F CARRIER: the carrier frequency of the target signal. 
The group speed used in calculating steering delays is 
based on this frequency. This must be entered in floating 
boant. 

@® FFT _ LENGTH: the length of the FFT used if frequency 
domain output is selected. This value must be a power of 
two which is less than or equal to the number of samples 
in a single channel of input data. 


@ SWAP: a macro used in calculating a FFT. 


The following program functions are adaptations of 


those found in Reference 11: 


@ polint: performs the polynomial interpolation. 
e realft: calculates Fast Fourier Transforms. 


@e fourl: calculates Fast Fourier Transforms 


oe 


@e vector: allocates memory for one dimensional arrays. 
@ matrix: allocates memory for two dimensional arrays. 
e free matrix: deallocates memory from two dimensional 
arrays. 
@ ExitOnError: provides abnormal program termination. 
2. Beamformer Source Code 
| kkk 
* PROGRAM: BEAMFORMER vsn 1.0 
* WRITTEN BY: Steven Crocker 
* LAST UPDATE: October 9, 1991 
* 
* This program takes input from various data files and the user. It 
x outputs a data file. The inputs are a number of channels of digital 
* acoustic data, and information regarding the physical characteristics 
* and geometry of the receiving array. Additionally, environmental data 
* in the form of normal mode eigenfunctions and eigenvalues at the 
* receiving array are required to operate this beamformer. The output 
* Ss a single channel of acoustic data. 
* 
takk / 
#define UNIX VERSION /* either ANSI or UNIX x/ 
Hdefine ASCII ON /* select output mode ON or OFF x / 
#define BINARY OFF /* select output mode ON or OFF x/ 
H#define SIGNAL ON /* either ON or OFF x/ 
H#define SPECTRUM ON /* either ON or OFF x/ 
H#define LOWER_SENSOR OFF /* either ON or OFF */ 
Hdefine VALIDATE OFF /* either ON or OFF x/ 
H#define ERROR_ESTIMATE OFF /* either ON or OFF */ 
#define ON 1 /* logical "switch" */ 
Hdefine OFF 0 /* logical "switch" */ 
H#define INTERPOLATE polint /* polint */ 
H#define ORDER 3 /* order of interpolator (odd) */ 
H#define STEP 1 /* number of steps for derivatives x/ 
H#define TINY 1.0e-25 /* prevents division by zero x/ 
#define PI 3.14159265359 /* for freq to omega conversions x/ 
H#define RADIAN 57.2957795131 /* for degree to radian conversions */ 
Hdefine OFFSET 4.0 /k dist btwn upper sensor and phone #1 */ 
Hdefine DELTA_R 45.0 /* array element spacing */ 
#define CTD_OFFSET 0.5 /*k diff btwn ctd depth inc & Ist depth */ 
H#define SSP_LENGTH 2500 /* max number of pts in eigunfunction %*/ 
H#define EIGVAL_LENGTH 230 /* max number of eigenvalues */ 
Hdefine LOOK_DIRECTION 217.0 /* direction from which signal arrives */ 
H#define TILT_BUFFER 120 /* max length of tilt data vectors */ 
H#define BUFFER TIME 12 /x input buffer length in seconds Cint)*/ 
H#define F_SAMPLE 228 /k sampling frequency Cint)*/ 
H#define CHANNELS 32 /* number of channels processed x/ 
H#define F_CARRIER 57.0 /* carrier frequency */ 
Hdef ine FFT_LENGTH 2048 /* radix 2 <= (BUFFER_TIME-2)*F_SAMPLE */ 
Hdefine SWAP(a,b) tempr=(a);(a)=(b); (b)=tempr 
#Hinclude<stdio.h> 
#include<malloc.h> 
#include<math.h> 
Hif defined ( ANSI ) 
#include<float.h> 
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#include<stdlib.h> 
int getInput (void); 
int putOutput (void); 
int processTilt(float **x, float **y, float **2z); 
int processModes(float **z, float *weight, float *ptrC); 
int dydx(float *x, float *y, float *ddx, int points); 
int polint(float *xa, float *ya, int n, float x, float *y, float *dy); 
int realft(float *data, int n, int isign); 
int fourl(float *data, int nn, int isign); 
float *vector(int length); 
float **matrix(Cint row, int col); 
int free_matrix(float **m, int row); 
void ExitOnError(char error_txtl]); 
Helif defined ( UNIX ) 
int vector(), 
matrix(); 
getInput(), ‘ 
putOutput(), 
processTilt(), 
processModes(), 
dydx(), 
polint(), 
realft(), 
fourl(), 
free_matrixQ), 
ExitOnError(); 
Hendif 


/k Global Variable Declarations */ 
float *xinSOUND, *outSOUND, *MeanSqError, Max=0.0, Min=1.0e25; 
int lastMinute, Minute=1, firstBuffer=1; 


main ) 
6 
FILE *fpVv1; 


amt 1, j),0K, Meeesmane. 


float *kx, **y, *kZ, *indx, *samples, arg, ans, err, xdelay, 
xdelta, *weight, *pwrSpectrum, Cgroup; 


/k Memory Allocation and Tilt Data Processing */ 
x=(float**)matrix(CHANNELS, TILT_BUFFER); 

y=(f loat**)matrix(CHANNELS, TILT BUFFER); 
z=(float**)matrixCCHANNELS, TILT BUFFER); 
processTilt(x, y, 2); 


if (VALIDATE==ON) 
{ 
printf("\nDumping array geometry to array.dat\n"); 


fpV1=fopen("array.dat","wt"); 
if (fpV1==NULL) ExitOnError("Error opening validation file"); 


TpeintiCipvi, Channet\t\t X\t\t Y\t\t Z\n"); 


for (1=1; 1<=lastMinute; i++) 
< 

fprintf(fpv1, "MINUTE: Z%i\n", 1); 

for (}=1; }<=CHANNELS; j++) 

SOR Inet CID veo tc. Af tet vce Ann. pe) ei) yj Ji), 20jICiI); 
} /* for */ 
fclose(fpv1); 
7k 1t */ 


/* Memory Allocationk/ 
weight=(float*) vector(CHANNELS); 
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indx=(f Loat*) vector (ORDER+1); 

delay=(f loat*) vector (CHANNELS) ; 

delta=(f loat*) vector (CHANNELS) ; 

out SOUND=(f loat*) vector ( (BUFFER TIME-1)*F_SAMPLE); 

INSOUND=(f Loat**)matrix(CHANNELS, BUFFER TIMEXF_ SAMPLE); 

if (Cshift=Cint*)mal loc( (CHANNELS+1)*sizeof (int) ) )==NULL) 
ExitOnError("Memory allocation failure for shiftC]."); 

for (1=1; 1<=ORDER+1; i++) indxliJ=(float)i; 

if CERROR_ESTIMATE==ON) 
MeanSqError=(f loat*) vector ( (BUFFER TIME-1)*F_SAMPLE); 
if (SPECTRUM==ON) pwrSpectrum=(f Loat*) vector (FFT_LENGTH/2) ; 

whi le(Minute<=LastMinute) 

{ 


/k Mode Data Processing */ 
processModes(z, weight, &Cgroup); 


/* Calculate delays, shifts, etc */ 

for (i=1; 1<=CHANNELS; 1++) 

{ 
delaylij=xfiJ£Minute]/Cgroup; /k Time delay = */ 
shiftlCij=Cint)(delaylij*x(float)F SAMPLE); /* # of samples */ 

/k fraction of 1 sample */ 

deltaLilJ=delayliJ*( float) F_SAMPLE-(float)shiftlidJ; 

} /* for */ 


if (VALIDATE==0N) 
fi 
printf("Dumping element delays to delay.dat\n"); 


if (firstBuffer) 


< 
if ((fpVv1=fopen("delay.dat", "wt")) == NULL) 
ExitOnError("Error opening validation file."); 
payee it */ 
else 


if ((fpVv1=fopen("delay.dat", "at")) == NULL) 
ExitOnError("Error opening validation file."); 
} /* else */ 
fprintf(fpV1,"MINUTE: %i\n", Minute); 


fprintt(fpv1, 
"Channel\t delay\t\t int shift\t fraction of 1 shift\n"); 


for (1=1; 1<=CHANNELS; i++) 


{ 
fprintt(fpv1,. "Zi\t Ae\te. ieee 
1,delayfij,shiftlCil,deltalid); 
+ /* for */ 
fclose(fpv1); 


printf("Dumping phone weights and grp speed to modal.dat\n"); 


if(firstBuffer) 
{ 
if ((fpV1=fopen("modal.dat", "wt")) == NULL) 
ExitOnError("Error opening validation file.\n"); 
fprintf(fpv1,"Group speed for F_CARRIER is: %g\n\n\n",Cgroup); 
} /k if */ 
else 
£ 
if ((fpV1=fopen("modal.dat","at")) == NULL) 
ExitOnError("Error opening validation file. \n"); 
fprintf(fpv1,"\nHydrophone weights for minute %i\n",Minute); 
} /* else if */ 


54 


fprintf(fpv1, "Channel \tWeight\n"); 


for (1=1; 1<=CHANNELS; i++) 
fprintf(fpv1, "“i\tZe\n",i1,weightlid); 


fclose(fpV1); 
Dark it */ 


1f (SPECTRUM==ON) 
for (k=1;k<=FFT_LENGTH;k++) pwrSpectrum(kJ=0.0; 


for (n=1; n<=60/(BUFFER_TIME-2); n++) 
i. 
getInput(); 


if(firstBuffer) /* Produce first output buffer */ 
{ 
for (i=1; i<=F_SAMPLE; i++) 
{ 
out SOUND(1)=0.0; 
if CERROR_ESTIMATE==ON) MeanSqErrorliJ=0.0; 
} /* for */ 
for (i=F_SAMPLE+1; i<=(BUFFER_TIME-1)*F_ SAMPLE; i++) 
{ 
outSOUND[1J=0.0; 
if CERROR_ESTIMATE==ON) MeanSqError[iJ]=0.0; 


for (j=1; j<=CHANNELS; j++) 
{ 

if (delay[j)]>=0.0) 

{ 


arg=(float) CORDER+1) /2.0+(1-deltalj]); 
samples = &inSOUND[jJCi-shift£jJ-CORDER+1) /2); 
Pale, it eR/ 
else if (delay[j1<0.0) 
{ 
arg=(f loat) CORDER+1)/2.0+deltalj]; 
samples = &inSOUNDCjJCi-shiftljJ-(CORDER+1) /2+1]; 
} (xtelse 1f */ 
INTERPOLATE(indx, samples, ORDER+1, arg, &ans, Serr); 
out SOUNDC1J=outSOUNDCiJj+ans; 
if (ERROR_ESTIMATE==ON) 
MeanSqErrorliJ=MeanSgErrorCil+errtxerr; 
} /* for */ 
if (fabs CoutSOUNDC1])>Max) Max=fabs(outSOUND[iJ); 
if (fabsCoutSOUND[LiJ)<Min) Min=fabs(outSOUND[i]); 
if CERROR_ESTIMATE==ON) 
MeanSqError[iJ=MeanSgErrorlil/(f loat) CHANNELS; 


} /* for */ 

» 7% it */ 

else /* Produce subsequent output buffers */ 

{ 
for (i=F_SAMPLE+1; i<=(BUFFER_TIME-1)*F SAMPLE; i++) 
{ 


outSOUND[i-F_SAMPLEJ=0.0; 
if (ERROR_ESTIMATE==ON) MeanSqError[i-F_SAMPLEJ=0.0; 


for (j=1; j<=CHANNELS; j++) 
{ 

if (delay[€jJ>=0.0) 

{ 


arg=(float) CORDER+1)/2.0+(1-deltaljJ); 
samples = &inSOUND[jIJCi-shiftljJ-CORDER+1) /2]; 
Dal Rot H/ 
else if (delay[jJ<0.6) 
{ 
arg=(float) (CORDER+1)/2.0+deltaljJ; 
samples = &inSOUND(j)JCi-shiftljJ-CORDER+1) /2+1]); 


a 


} /* else if */ 
INTERPOLATE(indx, samples, ORDER+1, arg, &ans, &err); 
out SOUNDLi-F_SAMPLEJ=out SOUNDL1-F_SAMPLEJ+ans; 
if (ERROR_ESTIMATE==ON) 
MeanSqErrorli-F_SAMPLEJ=MeanSqErrorli-F_SAMPLEJ+errxerr; 
> {® Tor */ 
if (fabs (out SOUND[ 1-F_SAMPLEJ)>Max) 
Max=fabs (outSOUNDLi-F_SAMPLE]); 


if (fabs (out SOUNDL1-F_SAMPLEJ) <Min) 
Min=fabs (out SOUND[1-F_SAMPLE]); 

if (ERROR_ESTIMATE==O0N) 
MeanSqErrorCli-F_SAMPLEJ=MeanSqError(li-F_SAMPLE]/ 

(float) CHANNELS; 
} /* for */ 
+} /* else */ 

if (SIGNAL==ON) putOutput(); 


if (SPECTRUM==ON) 

{ 
window(outSOUND, FFT_LENGTH); 
realft(outSOUND, FFT LENGTH/2, 1); 
for (k=0;k<FFT_LENGTH;k += 2) 
i 

pwrSpectrum[k/2+1J=pwr Spectrum[k/2+13+ 

out SOUNDLk J*xout SOUNDE kJ +out SOUNDEk+1]*outSOUND[kK+1]; 
+} /* for */ 

icf Neate 

firstBuffer=0; /* set firstBuffer to false */ 

+} /* for */ 


if (SPECTRUM==0N) 
{ 
for (k=1; k<=FFT_LENGTH; k++) 
pwrSpectrumLkJ=pwrSpectrum[kIJ* (float) (BUFFER_TIME-2) /60.0; 


dumpSpectrum(pwrSpectrum) ; 
} /* Vt 4/ 
printf("\t%1 minutes of input data processed.\n", Minute); 
Minutet+; /k increment minute counter */ 

+} /* while */ 
printfC"EXECUTION COMPLETE: End of tilt data encountered\n"); 
printf('"Maximum magnitude encountered was: Z%e\n",Max); 
printf("Minumum magnitude encountered was: %e\n",Min); 
exit(0); 

FL RAKIIK KIA IKK AIK IAI IA IAI IIIS KI KIKI IIIS ISA SISA IIIA IIIS ASAI AAA AAAI 
DI III END oman odode / 
[RKKKKKEKKKISAK KEKE KI KER KKK EEK IK IAI RRR RAIS AISI IAI AAAI RII AIIAAIA / 

[ xkKkr 


* FUNCTION: getInput() 

x 

x This function handles all acoustic input. It also provides one of 
x two normal process terminations available in the program. (The other 
* is located in main().) 

x 

* Arguments: none 

* 

* Return value: 0 

x 

* Functions called: vector () ExitOnError() 

x 

* Definitions called: ANSI UNIX 

* F_SAMPLE CHANNELS 

* BUFFER TIME 

x 

* Global variables called: inSOUNDTICI Min 

* firstBuffer Max 

x 
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x 
x 


Significant memory allocation: diskBufferl] 


wrk / 

#if defined ( ANSI ) 
int getInput (void) 
#Helif defined ( UNIX ) 
get Input () 

Hendif 


t 


int i, j, buffer, items; 
float *diskBuffer; 

char fileName[80]; 
static FILE *fpStatic; 


if (firstBuffer) 
tt 
printf("Enter file name for input acoustic data: "); 


if((scanf("%s", fileName) )==EOF) 
ExitOnError("Fatal error in scanf()"); 


Brame \n\n: 2; 


if ((fpStatic=fopen(fileName, "rb")) == NULL) 
ExitOnError("Error opening INPUT ACOUSTIC data file."); 
eye it x/ 


/x* Memory Allocation */ 

if(firstBuffer) buffer=BUFFER TIME*F_SAMPLEXCHANNELS; 
else buffer=(BUFFER TIME-2)*F SAMPLEXCHANNELS; 
diskBuffer=(float*) vector (buffer) ; 


items=fread((char*) (diskBuffer+1), sizeof(float), buffer, fpStatic); 

if (items==buffer) /xcontinue*/; 

else if(ferror(fpStatic) != 0) 
ExitOnError("Error encountered while reading input acoustic data"); 

else if(feof(fpStatic) != 0) 

{ 
DPI Nth CNAN\ TERRIER AHA RK AKAIKE III IIIA KR IIA IIIA IR ISA IIS IAAI IIIB PY 
Drintt ('\t End of File reached: EXECUTION COMPLETE\n"); 
printf("\t\tZ1 minutes of data processed\n",Minute-1); 
printf("\t\t41 bytes of data discarded\n", itemsxsizeof(float)); 
printf("\t\tMaximum magnitude encountered was: %e\n",Max); 
printf("\t\tMinimum magnitude encountered was: Z%e\n",Min); 
Brantt("\t End of File reached: EXECUTION COMPLETE\n"); 
DRINT ECU \ TERRA KKKKAAHK ARKH KK KIKKKK KKK AAA IAEEAAARIARARAAIIAE EK \ 1) 2 
fclose(fpStatic); 
exit(0); 

} /* else if */ 

else ExitOnError("Unknown error handling acoustic input file."); 


for (i=1; 1<=BUFFER_TIMEXF SAMPLE; i++) 
{ 

for (j=1; j<=CHANNELS; j++) 

de 


if (firstBuffer) 


inSOUNDEj ICI = diskBuf ferl[CHANNELS*(i-1)+j]; 
Fayk it */ 
else 
- 
if (1<=2*F_ SAMPLE) 
INSOUNDE jI£II=iMSOUNDL)ICi+(BUFFER_TIME-2)*F_SAMPLEJ; 
else 
InSOUNDE jICiJ=diskBuf fer [CHANNELS*(i-2*F_SAMPLE-1)+j]; 
> /* else x/ 
> /* for */ 
ya * tor */ 
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/* Deallocate Memory */ 
free((char*)diskBuffer); 
return 0); 

> aaa / 
/kkxekek END getInput **x%*/ 
Doo oo oak / 


tI / 
Hif defined ( ANSI ) 
Int putOutput (void) 
Helif defined ( UNIX ) 
putOutput () 
Hendif 
di 
mt 15 eur 
char fileName[12], mode[2]; 


static FILE *fpOutSound, *fpMSE; 


if (firstBuf fer) 
{ 
if (ASCII==ON) 
{ 
mode([OJ='w'; 
mode([1J='t'; 
alk it 87 
else if (BINARY==ON) 
4 
mode(OJ='w'; 
mode(1J='b'; 


> /k else if */ 


cut=1; 


[kkekk 

* FUNCTION: putOutput() 

* 

* This function handles all acoustic output. Additionally, it outputs 
x the estimated mean squared error from the interpolators (if enabled). 
* 

* Arguments: none 

* 

* Return value: 0 

* 

* Functions called: ExitOnError() 

* 

x Definitions called: ANSI UNIX 

x F_SAMPLE BUFFER TIME 

* 

* Global variables called: out SOUND) MeanSqError(] 

* firstBuf fer 

x 

* Significant memory allocation: none 

x 


printf("Enter file name for output data: "); 


if((scanf("%s", fileName) )==EOF) 
ExitOnError("Fatal error in scanf()"); 


printf("\n\n"); 


if ((fpOutSound=fopen(fileName, mode)) == NULL) 
ExitOnError("Error opening OUTPUT data file."); 


if (ERROR _ESTIMATE==ON) 


i 
printf("Opening file error.dat\n\n"); 
if (CfpMSE=fopen("error.dat", mode)) == NULL) 
ExitOnError("Error opening error.dat"); 
} /* if */ 
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Pek it */ 
else cut=2; 


if CERROR_ESTIMATE==ON ) 
< 

if CASCII==0N) 

‘ 


for (i=1; 1<=(BUFFER_TIME-cut)*F_SAMPLE; 1++) 
fprintf(fpMSE,"%e\n", MeanSgErrorlil); 
Beye if */ 
else if (BINARY==ON) 
< 
if (fwrite((char*) (MeanSqError+1),sizeof (float), (BUFFER _TIME-cut)* 
F_SAMPLE, fpMSE)==(unsigned) (BUFFER_TIME-cut)*F_SAMPLE) ; 
else if(ferror(fpMSE) '= Q) 
ExitOnError("Error encountered writing error data"); 
else ExitOnError("Unknown error handling error file."); 
>} /k else if */ 


} /k if */ 
if (ASCII==0N) 
{ 
for (1=1;1<=(BUFFER_TIME-cut)*F_SAMPLE; i++) 
fprintf(fpOutSound, "Z%e\n", outSOUND[iJ); 
cere it */ 
else if (BINARY==ON) 
e 


if (fwrite((chark) CoutSOUND+1),sizeof (float), (BUFFER TIME-cut)* 
F_SAMPLE, fpOutSound)==(unsigned) (BUFFER _TIME-cut)*F_SAMPLE) ; 
else if(ferror(fpOutSound) != 0) 
ExitOnError("Error encountered writing output acoustic data"); 
else ExitOnError("Unknown error handling acoustic output file."); 
} /* else if */ 
return( QO ); 
FL KKIAKKKK KR KKK RK KKKEERE RE / 
[kkkkKK END putOutput ***xe/ 
[KKKIKAKEKKIEKKRIAKIA IAAI / 
[xkkkk 


A FUNCTION: dumpSpectrum() 

: This function handles dumps the signal power spectrum (if selected). 
* Arguments: pwrSpectrum 

Return value: 0 

: Functions called: ExitOnError() 

* Definitions called: ANSI UNIX 

* F_SAMPLE FFT_LENGTH 

; Global variables called: none 

; Significant memory allocation: none 

Pe, 


Hif defined ( ANSI ) 
int dumpSpectrum ( float *pwrSpectrum ) 
Helif defined ( UNIX ) 
dumpSpectrum ( pwrSpectrum ) 
float *pwrSpectrum; 
Hendif 
£ 
vine: 
static float sequence=1.0; 
char fileName[12], mode[2]; 
static FILE *fp; 


a0 


if (firstBuffer) 


t 
if (ASCII==0N) 
{ 
mode[OJ='w'; 
mode(1]='t'; 
>} /k if */ 
else if (BINARY==ON) 
fi 
mode{OJ='w'; 
mode{1J='b'; 
} /*k else if */ 
printf("Enter file name for output SPECTRUM data: "); 
if((scanf("%s", fileName) )==EOF) 
ExitOnError("Fatal error in scanf()"); 
Printt. \nvne): 
if ((fp=fopen(fileName, mode)) == NULL) 
ExitOnError("Error opening OUTPUT SPECTRUM data file."); 
printf("Select output format: \n"); 
printf(" Enter O for MATLAB compatible output. \n"); 
printf(" Enter 1 for GRAFTOOL compatible output.\n"); 
if((scanf("%1", &1))==EOF) 
ExitOnError("Fatal error in scanf()"); 
if(¢i '= 0) && Ci $= 1)) ExitOnError("Invalid output selection"); 
>} /k if */ 
if (ASCII==0N) 
{ 
if(firstBuffer && format==1) 
< 
forintt (ip, 0. ~); 
for (1=1, 1<=FFT_LENGTH/2; i++) 
fprintf(fp,"%e ", (float) (i-1)*F_SAMPLE/FFT_LENGTH) ; 
Tpeintt ip; me: 
Peet X77 
if (format==1) fprintf(fp,"%e ",sequence); 
for (1=1;1<=FFT_LENGTH/2; i++) 
fprintf(fp,"%e ", pwrSpectrum[iJ); 
ToORINtICte,. \Me 
sequence=sequencet1.0; 
Pex eat x) 
else if (BINARY==ON) 
{ 


if (fwrite((charx) (pwrSpectrum+1),sizeof (float),FFT_LENGTH/2, fp)== 
(unsigned) FFT_LENGTH/2) ; 
else if(ferror(fpOutSound) != 0) 
ExitOnError("Error encountered writing output SPECTRUM data"); 
else ExitOnError("Unknown error handling acoustic SPECTRUM file."); 
>} /k else if */ 
return( O ); 
3 add od kkk / 
[xkkkk end dumpSpectrum *xxrk/ 
[oI In / 


| kxkke 

* FUNCTION: processTilt() 

* 

* This function handles all array tilt data. It calculates the X, Y, Z 
x coordinates of each hydrophone as a function of time. The coordinate 
* system is oriented such that X points toward the signal "origin" and 
* 2Z points down. 

* 
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% >t Ot OF OF OF OF OF OF OF ot OF ot oF Ot tot ot ot Ot OF 


Arguments: xCJCI yCjC) 


ZEAba 
Return value: 0 
Functions called: ExitOnError() matrix() 
vector () free_matrix() 
Definitions called: ANSI UNIX 
DELTA_R LOWER SENSOR 
RADIAN CHANNELS 


LOOK DIRECTION OFFSET 
TILT_BUFFER 


Global variables called: lastMinute 
Significant memory allocation:  xx(CICJ yyCIJCI 
ZZE DE) Ciel) 
anglel] udepth(] 
ldepthl] 
tok / 


Hif defined ( ANSI ) 

int processTilt(float **x, float **y, float **z) 
Helif defined ( UNIX ) 

processTilt(x,y,2Z) 

float *kx, Key, **zZ; 

#endif 


i‘ 


ime), 3, NotEOF; 

float *xtilt, *angle, *kudepth, *ldepth, **xx, *kyy, **zz, theta; 
char fileNamel12]; 

FILE xfp1, *fp2; 


/* Open Data Files */ 
printf("Enter file name for upper tilt sensor data: "); 


faccscant("“s", fileName) )==EOF) 
ExitOnError("Fatal error in scanf()"); 


printf("\n\n"); 


if ((fpl=fopen(fileName, "rt")) == NULL) 
ExitOnError("Error opening UPPER TILT data file."); 


if (LOWER_SENSOR==ON) 
< 


printf("Enter file name for lower tilt sensor data: "); 


if((scanf("%s", fileName) )==EOF) 
ExitOnError("Fatal error in scanf()"); 


printt C"\n\ni); 


if (Cfp2=fopen(fileName, "rt")) == NULL) 
ExitOnError("Error opening LOWER TILT data file."); 
idepth=(float*) vector (TILT_BUFFER) ; 
Poet Xf 


/x Memory Allocation */ 
tilt=(float*)vector(TILT_BUFFER); 
angle=(float*) vector(TILT_BUFFER) ; 
udepth=(float*) vector (TILT_BUFFER) ; 

xx=(f loat**)matrix(CHANNELS, TILT BUFFER); 
yy=(f loatkk)matrix(CHANNELS, TILT BUFFER); 
zz=(f loatk*)matrix(CHANNELS, TILT BUFFER); 


i=1; /k read upper tilt data */ 
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notEOF=1; 
whileCnotEOF ) 
{ 


if (fscanf(fp1,"%g %g “%g\n",Stiltlil,&angleli),&udepthlid) != EOF) i++; 
else notEOF=0; 


} 
if (LOWER_SENSOR==ON) 
£ 
j=), /* read lower tilt data */ 
NOtEOF=1; 
whi leCnotEOF ) 
x 
if (fscanf(fp2, "%g\n",&ldepthl[jJ>) != EOF) j++; 
else notEOF=0; 
+ /* while */ 
if (i<=j) lastMinute=1-1; 
else lastMinute=j-1; 
7 


else lastMinute=i-1; 


[&KKKKKKKKKKKTHIS 1s the assumed array geometry: LINEARKKKkkkkkkKx / 
for (j=1; j<=CHANNELS; j++) 
{ 


for (1=1; 1<=lastMinute; 1++) 
{ | 
xxCjJDIJ=DELTA_R*( float) (j-1)*cos(angleliJ/RADIAN)* 
sin(tiltliJ/RADIAN) ; 
yyljICIJ=DELTA_R*(f Loat)(j-1)*sinCangleliJ/RADIAN)* 
sin(tilt£1J/RADIAN) ; 
ZZ20jJIJCIJ=DELTA_R*( float) (j-1)*cos(tiltlil /RADIAN)+ 
OFFSET*cos(tiltliJ/RADIAN)+tudepthliJ; 
} /* for */ 
y 1% for x7 
LEKKKKKKIIIKKAAA III RIA IKARIA IATA KERRIER AAAAAAI KE / 


theta=(360.0-LOOK_DIRECTION)/RADIAN; /* coordinate rotation */ 
for (j=1; j<=CHANNELS; j++) 
< 


for (1=1; i<=lastMinute; 1++) /k points x into signal */ 
£ 
xCjJCiJ=xxCjICil*cos(theta)-yy[jICil*sin(theta); 
yCjJCij=xxCjllilxsin(thetad+yyljJCil*cos (theta); 
z(jJCiJ=zz(jICil; 
+ /* for */ 
> /* tor x7 


/x Memory Deallocation */ 
if (LOWER_SENSOR==ON) free((charx)ldepth); 
free((char*®) tilt); 


free((char*)angle); 
free((char*)udepth); 
free_matrix(xx, CHANNELS); 
free_matrix(yy, CHANNELS); 
free_matrix(zz, CHANNELS); 
fclose(fp1); ~ 
if CLOWER_SENSOR==ON) fclose(fp2); 
return( 0 )>; 
3 SI od / 
/kkkkk END processTilt *kii*/ 
id oo noo i / 
[kkekk 
* FUNCTION: processModes( ) 


This function handles the normal mode data. It calculates hydrophone 
weights and group speed. The user must insure that the depth vector 
and eigenfunction vector are of equal Length. 


+ 4+ 4+ 
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* 

* Arguments: Zt weightl) 
x ptrc 

x 

* Return value: 0 

x 

* Functions called: vector () ExitOnError() 
k INTERPOLATE () dydx() 

* 

x Definitions called: ANSI UNIX 

* PI SSP_LENGTH 
x EIGVAL_LENGTH ORDER 

x CHANNELS F_CARRIER 
* 

x Global variables called: Minute 

x 

* Significant memory allocation: depth] Zn) 

x Oonofflj wo 

x Kr{J dwdK 

* 

kik / 


Hif defined ( ANSI ) 

int processModes(float **z, float *weight, float *ptrc) 

Helif defined ( UNIX ) 

processModes(z,weight,ptrc) 

float **2z, kweight, *ptrC; 

Hendif 

£ 
int i, j, ptsEigVal, set, notEOF, deadPhones, weightNotAssigned; 
static int ptsEigfun; 
float %w, *Kr, *dwdK, err; 
static float depth[SSP_LENGTH+1], Zn[€SSP_LENGTH+1], OnOffCCHANNELS+1]; 
char key, fileNamel12]; 
FILE *fp1, *fp2; 


if (firstBuffer==1) 

< 
w=(f loat*) vector(EIGVAL_ LENGTH); 
Kr=(f loat*) vector (EIGVAL_LENGTH) ; 
dwdk=( f loat*)vector(EIGVAL_LENGTH) ; 


printf("Enter file name for normal mode data (eigenfunction): "); 


if((scanf("%s", fileName) )==EOF) 
ExitOnError("Fatal error in scanf()"); 


print fC"\n\n"); 


if ((fpl=fopen(fileName, "rt")) == NULL) 
ExitOnError("Error opening NORMAL MODE data file (eigenfunction)"); 


printf(C"Enter file name for normal mode data (eigenvalues): "); 


if((scanf("%s", fileName) )==EOF) 
ExitOnError("Fatal error in scanf()"); 


prints C\n\n"); 


if ((Cfp2=fopen(fileName, "rt")) == NULL) 
ExitOnError("Error opening NORMAL MODE data file (eigenvalues)."); 


1=1; /® read normal mode data */ 

notEOF=1; 

while (notEOF ) 

a 
if(fscanf(fpi, "%g %g \n", &depthlil, &ZnCiJ) '= EOF) i++; 
else notEOF=0; 

- 
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ptsE1igFun=i-1; 
for(i=1; 1<=ptsEigFun; i++) depth£Cil=depth£ij+CTD_OFFSET; 


1=1; 

notEOF=1; 

while (notEOF) 

+ 
if(fscanf(fp2, "%g %g \n", &wlil, &KrOij) != EOF) i++; 
else notEOF=0; 

} 

ptsEigVal=i-1; 


for (i=1; 1<=CHANNELS; i++) OnOffCiJ=1.0; 
printf("Do you want to turn off any hydrophones? "); 


if ((scanf("Z%s",&key) )==EOF ) 
ExitOnError("Fatal error in scanf()"); 


if(key=="y' ,, key=='Y') 
< 
printf("\nHow many hydrophones must be secured? "); 


if((scanf("Z%1", &deadPhones) )==E0F ) 
ExitOnError("Fatal error in scanf()"); 


for(1=1; 1<=deadPhones; i++) 
t 


printf("\nEnter hydrophone number to secure: "); 


if((scanf("%Zi", &j))==EOF) 
ExitOnError("Fatal error in scanf()"); 


if (j>CHANNELS |; j<1) 
ExitOnError("Bad hydrophone identification"); 
Onof f0jI=0.0; 
} /* for */ 
} /k if */ 
[Rts / 


for(1=1; i<=lastMinute; i++) 
t 
if (depth[pt sEigFunJ<z{CHANNELSIJCiJ) 
{ 
printf("Max eigenfunction depth is: %f\n",depth[ptsE1gFun]); 
printf("at depth[] index of: %i\n",ptsEigFun); 
printf("Max depth of phone number %1 1s: %f\n\n", 
CHANNELS, zZDCHANNELSJ(1); 
ExitOnError("Fatal data set error"); 
} /k if */ 
} /* for */ 


jay 

for (€1=1; 1<=CHANNELS; i++) 

< 
we ightNotAssigned=1; 
while (j<=ptsEigFun && weightNotAssigned) 
< 


if (zDiJ[Minutel<0.0 |; depthlj1<0.0) 
{ 
print tC i=Zz1\n", 12; 
printf("j=Z1\n",)); 
printf (C"Minute=Zi1\n" Minute) ; 
printf("zCiljCMinute] is: %1r\n", zCil{Minute)); 
printf("depthlj] is: “f\n",depth{jJ); 
printf("Depth less than zero encountered in processModes."); 
Print Ca \a\m. 2; 
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} 


ExitOnError("Check input depths for coordinate orientation"); 
D/% 1f */ 


if(depth(jJ<zCiJ(Minute] && depthlj+1]>zCiICMinutel) 
{ 


set=(ORDER+1)/2; 
INTERPOLATE(&depth[j-setJ, &Zn[Cj-set], ORDER+1, zCiJ(Minutel, 
&weightliJ, &err); 

weightNotAss igned=0; 

Dee it */ 

let, 

} /* while */ 
} /* for */ 


for (i=1; i<=CHANNELS; i++) weightliJ=On0fflij*weightliJ; 


if (Minute==1) 
{ 
dydx(Kr,w,dwdkK,ptsEigVal); 
for (i=1;1<=ptsEigVal; i++) 
c 
if (wlLid<2.O*%PI*F_CARRIER && wlit+1]>2.O*%PI*F_CARRIER) 
{ 
set=(ORDER+1)/2; 
INTERPOLATE (&wli-set],&dwdkCi-set],ORDER+1, 
2.Q*%PIXF CARRIER, ptrC,&err); 
} /k if */ 
} /* for */ 


free((chark)w); 
free((char*)Kr); 
free( (char )dwdk); 
> /k if */ 
return( 0 ); 
[SAKIAAIKKEKAKIAAIIIATIAAKIKRE / 


/kxk «kk END processModes ***k*/ 
[KK KIA KEIK KEIR IKI AAA / 


[keer 

* FUNCTION: dydx() 

* 

* This function estimates derivatives. 

* 

* Arguments: xC) yCJ 
* ddxlJ points 
* 

* Return value: 0 

* ; 

* Functions called: ExitOnError() 

x 

* Definitions called: ANSI UNIX 
x STEP 

* 

* Global variables called: none 

* 

* Significant memory allocation: none 

* 

KKK / 


#if defined ( ANSI ) 

int dydx(float *x, float *y, float *ddx, int points) 
Helif defined ( UNIX ) 

dydx(x,y,ddx, points) 

float *x, *y, *ddx; 

int points; 

fendif 


{ 


nt 1; 
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for (n=1;n<=points;n++) 
{ 
if ((n>=STEP) && (n<=points-STEP)) /*center*/ 
ddxEnJ=(y(n+STEPJ-y[n-STEP])/(xEn+STEPJ-xEn-STEP J); 


else if (n<STEP) /*beginning*/ 
ddx(nJ=(y(n+STEPJ-y£1})/(xEn+STEPJ-x01]); 


else if (n>points-STEP) /kend*/ 
ddx(nj=(yCpointsj-y£n-STEP])/(xLpointsJ-xEn-STEPI); 


else 
ExitOnError("Index error in dydx"); /* sanity check */ 
} /* for */ 
return( O ); 
3 Sionoinobik / 
[kkk END dydx *kikk/ 


Dikhit / 
| xkickk 
* FUNCTION: polint() 
x 
* This function performs polynomial interpolation. 
* 
* Arguments: xaCl] yal) 
x n x 
x y dy 
x 
* Return value: 0 
x 
* Functions called: vector () ExitOnError() 
r 
x Definitions called: ANSI UNIX 
x 
x Global variables called: none 
* 
* Significant memory allocation: df] ety 
x 
KKKKK/ 


Hif defined ( ANSI ) 
int polint( float *xa, float *ya, int n, float x, float *y, float *dy) 
Helif defined ( UNIX ) 
polint(xa,ya,n,x,y,dy) 
float *xa, *ya, x, ky, *dy; 
inten: 
Hendif 
{ 
Int Wem, ns=1- 
float den, dif, dift, ho, hp, w; 
float *c, *d; 


dif=fabs(x-xal[1]); 


c=(f loat*)vector(n); 
d=(float*)vector(n); 


for (1=1; 1<=n; 1++) 


if ((dift=fabs(x-xa[i])) < dif) 
{ 
NS=1; 
dif=dift; 
} /k if */ 
cCij=yalid; 
dCij=yalil; 


} /* for */ 
ky=yalns--]; 

for (m=1; m<n; m++) 
x 
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for (1=1; 1<=n-m; i++) 
- 
ho=xaflij-x; 
hp=xali+m]-x; 
w=cCi+1)-dCiJ; 
if ((den=ho-hp)==0.0) 
ExitOnError("Error in routine POLINT"); 
den=w/den; 
d({iJ=hp*den; 
clij=hox%den; 
+ /* for */ 
ky += (xdy=(2kns < (n-m) ? c£ns+1] : dl€ns--])); 
} /* for */ 
free((char*)d); 
free((chark)c); 
return( 0 ); 
3 Dd ook / 
[kkkkK END polint ****/ 


[IIIT TAA TIAA / 

[keekk 

* FUNCTION: window() 

* 

* This function applies a Blackman window to a vector. 
x 

* Arguments: dataC] N 
* 

* Return value: 0 

x 

* Functions called: none 

x 

x Definitions called: ANSI UNIX 
* PI 

x 

x Global variables called: none 

x 

* Significant memory allocation: none 

x 

kkeek / 


#if defined ( ANSI ) 

int window(float *data, int N); 
Helif defined ( UNIX ) 

window( data, N ) 

float *data; 

int N; 

Hendif 

if 


Tat; 


for (n=0; n<N; n++) 
< 
data[n+1J=data[n+1J]*(0.42+0.5xcos (2. O0*xPI* (float) (n-N/2)/(f Loat) (N-1)) 
+0.08%cos(4.O*PI*( float) (n-N/2)/( float) (N-1))); 
} /* for */ 
return ( 0 ); 
3 Dando / 
[kkkkk end window **rk* / 
[KAAKIAAAIAAAAAAAA AAAI / 


[keekk 

* FUNCTION: realft() 

* 

* This function calculates FFT's 

* 

* Arguments: dataC] n 
* isign 

* 

* Return value: 0 

* 
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zs Functions called: four1() 

x Definitions called: ANSI UNIX 
Global variables called: none 

* Significant memory allocation: none 

rei 


Hif defined ( ANSI ) 

int realft(float *data, int n, int isign) 
Helif defined ( UNIX ) 

realft(data, n, isign) 

float *data; 


Int 
Hendif 
{ 


isign; 


Int (1), 117.278 1555147 eps, 
float cl=0:5, c2, hir, hli, ner, hei; 
double wr, wi, wpr, wpi, wtemp, theta; 


theta=3.141592653589793/ (double)n; 
if Cisign==1) 


1 
ce = -0.5; 
fourl(data, n, 1); 
> /e if */ 
else 
{ 
c2=0.5; 


theta = -theta; 
> /* else */ 
wtemp=sin(0.5xtheta) ; 
wpr = -2.0*wtemp*wtemp; 
wpi=sin(theta) ; 
wr=1.Q+wpr; 
Wi=wpi; 
n2p3=2kn+3; 
for (122; i<=n/2; i++) 


{ 
14=1+ (13=n2p3-(12=14 (11=14+1-1))); 
hir=ci*(datalilj+datali3)); 
h1i=ci*(datali2l]-datali4)); 
her = -c2x(datali2j+datali4]); 
h2i=c2ek(dataliij-datali3]); 
datalilj=h1r+we kher-wikh2i; 
datali2J=h1itwrth2itwikher ; 
datali3J=h1r-wrther+wikh2i ; 
datali4] = -hlitwrth2itwikher; 
wr=(wtemp=wr ) kwpr-wi*wpitwr; 
Wi=wikwpr+wtemp*wpi+wi; 

> /* for */ 

if (isign==1) 

a 
data(1J=(hir=datal1))+datal2]; 
data(2)=hir-datal2]; 

> /k if */ 

else 

{[ 


datal1J=c1*( (hi r=datal13)+datal2]); 
datal2J=c1*(hir-datal2)); 
fourl(data,n,-1); 

} /k else */ 

return ( QO ); 


FL RAK KIKI KIKI III / 
[kkkxkk end realft *xtKx/ 
[KKIKKIKKIIKKK AREER IAI / 


|[ xkkKrk 
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* FUNCTION: four1() 

x 

* This function calculates FFT's 

x 

* Arguments: dataCl nn 
x isign 

x 

* Return value: 0 

x 

* Functions called: none 

x 

* Definitions called: ANSI UNIX 
x SWAP 

x 

* Global variables called: none 

x 

* Significant memory allocation: none 

x 


KkKK / 

#if defined ( ANSI ) 

int fourl(float *data, int nn, int isign) 

Helif defined ( UNIX ) 

fourl(data, nn, isign) 

float *data; 

int nn, isign; 

#endif 

£ 
int n, mmax, m, J, istep, 1; 
double wtemp, wr, wpr, wpi, wi, theta; 
float tempr, tempi; 


n=nn << 1; 
j=1; 
for (1=1; 1<n; i+=2) 
a 
lveeiee 1) 
< 
SWAP (dataljJ,datalil); 
SWAP (datal j+1),datalit+1)); 
>} /* for x/ 
m=n >> 1; 
while (m >= 2 && j > m 
{ 
iar 0; 
m >>=1; 
} /*® while */ 
tet 
+} /* for */ 
mmax=2; 
while (n > mmax) 
< 
istep=2%mmax; 


theta=6.28318530717959/(isignkmmax) ; 
wtemp=sin(0.5*theta); 
wpr = -2.0*wtemp*wtemp; 
wpi=sin(theta); 
wr=1.0; 
wi=0.0; 
for (m=1; m<mmax; m+=2) 
{ 
for (i=m; i<=n; i+=istep) 
1 
j=1+mmax ; 
tempr=wrxdataljJ-wixdatalj+1]; 
tempi=wrkdatalj+1]+witdatalj J; 
datal[jJ=datalij-tempr; 
data{j+1J=datalit+t1)]-tempi; 
datafi] += tempr; 
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dataCi+1) += tempi; 
ae One x7 
wr=(wtemp=wr )*kwpr-wikwpitwr; 
Wi=wixwpr+twtempxwpitwi; 
y 7* ton k/ 
mmax=1step; 
} /k while */ 
return (0 ); 
3 Dodi: / 
[xkkkk end four) *kxkk/ 
LRKKKKKKKKEKRAEKKAAEIEE / 


[kkk 

5 FUNCTION: vector() 

3 This function allocates memory for UNIT OFFSET vectors. 
* Arguments: Length 

Return value: ky 

* Functions called: ExitOnError() 

: Definitions called: ANSI UNIX 
; Global variables called: none 

Significant memory allocation:  vfC] 


kkk / 
Hif defined ( ANSI ) 
float *vector(int length) 
Helif defined ( UNIX ) 
vector( length) 
int length; 
Hendif 
‘ 

float xv; 


if ((v=(float*)malloc((length+1)*sizeof (float) ))==NULL) 
ExitOnError("Memory allocation failure in vector()."); 

#if defined ( ANSI ) 

return v; 
Helif defined ( UNIX ) 

return (long int)v; 
H#Hendif 
3 AI AIK KI / 

/xkxekk END vector kxkkk/ 

[SAKA AKIAIAI IIA IIIA AA / 


[kkkKk 

* FUNCTION: matrix() 

; This function allocates memory for UNIT OFFSET 2-D arrays. 
; Arguments: row col 

: Return value: km 

: Functions called: ExitOnError() 

: Definitions called: ANSI UNIX 

: Global variables called: none 

: Significant memory allocation: m(JC] 

x 


ikke / 


Hif defined ( ANSI ) 
float **xmatrixCint row, int col) 
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felif defined ( UNIX ) 
matrix(row,col) 
imt row, col; 
#endif 
{ 
int 1; 
float **m; 


if ((m=(float**) malloc ( (unsigned) Crowt+1)*sizeof (float*) ))==NULL) 
ExitOnError("Allocation failure 1 in matrix ()"); 
for (i=1; i<=row; i++) 
a 
if ((mCiJ=(float*)malloc( (unsigned) (col+1)*sizeof (float) ) )==NULL) 
ExitOnError("Allocation failure 2 in matrix()"); 
>} /* for */ 
#if defined ( ANSI ) 
return m™, 
Helif defined ( UNIX ) 
return (long int)m; 
Hendif 
Dd LRKKKKKKKKEEEEEKEEE KKK / 
[/kkkkx END matrix *xxxx/ 
[KKKKEKKAAKIAKKEKIK AEE / 
[eek 


FUNCTION: free_matrix( 

* This function deallocates memory from UNIT OFFSET 2-D arrays. 
; Arguments: mC ICI row 

‘ Return value: 0 

: Functions called: none 

; Definitions called: ANSI UNIX 

: Global variables called: none 

* 

* Significant memory allocation: none 

ney, 


Hif defined ( ANSI ) 

void free_matrix(float **m, int row) 
Helif defined ( UNIX ) 
free_matrix(m, row) 

float **m; 

Int row; 

Hendif 

i 


anit. 1; 


for(i=row; 1>=1; i--) 
free((char*)mliJ); 
free((charx)m); 
return( OQ ); 
3 Dag dick iok / 
/kkkkk END free matrix *kkk*/ 
LKKKKKKEKKKKRKKEEKKEEKKEEEKK IKE / 


[kkk 

* FUNCTION: ExitOnError() 

* 

* This function performs an abnormal process termination. 
* 

* Arguments: error_txtl] 

* 

i Return value: none 

* Functions called: none 
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; Definitions called: ANSI UNIX 
Global variables called: none 

* Significant memory allocation: none 

ee 


Hif defined ( ANSI ) 

void ExitOnError(char error_txtl]) 
Helif defined ( UNIX ) 
ExitOnError(error_txt) 

char error_txt[]; 


#Hendif 

{ 
fprintf(stderr,"Program run-time error ...\n"); 
fprintf(stderr,"%s\n",error_txt); 
fprintf(stderr,"...now exiting to system...\n"); 
exit(0); 


SL RRKKKKKAAAAAKKKEKKEEEEEKEKE / 


[kkkkR END ExitOnError «kerk / 
[LRKAKRAKIAAAIAIA IAAI AIK / 


B. DECOMMA.C 
This program removes commas and colons from the output of 


the upper tilt sensor of the Heard Island Array. 


[kkk 

* PROGRAM DECOMMA.C vsn 1.0 

* WRITTEN BY: Steven Crocker 

* LAST UPDATE: August 3, 1991 

x 

x This program takes any ASCII file and copies it to a user 

x defined file. It removes the commas and colons. 

Wak KK / 

#define VERSION ANSI /x either ANSI or UNIX x/ 
#define C 25 /k max Length for file names */ 
H#define COMMA 44 /* ASCII codes xf 


#define COLON 58 
Hdefine SPACE 32 
Hinclude<stdio.h> 
Hinclude <malloc.h> 


#if (CVERSION==ANSI) 
#include<stdlib.h> 

void ExitOnError(char error_txt[]); 
endif 


main() 
{ 
FILE *sacm_fp, *out_fp; 
char c, sacm_file[CJ, out_file([C]; 


printf("\n\n\nEnter file name for input SACM data. \n"); 
scanf("%s",sacm_file); 
Or imt tC Nana 2, 
printf("Enter file name for output SACM data.\n"); 
scanf("%s",out_file); 
if ((sacm_fp=fopen(sacm_file, "rt")) '= NULL) { 
if (Cout_fp=fopen(out_file, "wt")) != NULL) { 
while((c=fgetc( sacm_fp )) != EOF) { 


eZ 


if((c==COMMA) |! (c==COLON)) c=SPACE; 
TPULCCE, Our Tp), 
> /*k while */ 
} /* if */ 
else ExitOnError("Error opening output data file"); 
} /k if */ 
else ExitOnError("Error opening input data file"); 
fclose(sacm fp); 
fclose(out_fp); 
exit(0); 
} 


Hif (VERSION==ANSI) 

void ExitOnError(char error_txt[J) 
Helif (VERSION==UNIX) 
ExitOnError(error_txt) 

char error_txtC]; 


Hendif 

{ 
fprintf(stderr,"Program run-time error ...\n"); 
fprintf(stderr, "%s\n",error_txt); 
fprintf(stderr,"...now exiting to system...\n'"); 
exit(1); 

} 


C. SACM1.C 

This program reads the data output by the upper tilt 
sensor and formats it for use in the beamformer. Le 
calculates 60 second averages for all data fields, converts 
pressure to depth and calculates tilt direction based on 


eoerent velocity. 


| kkkkk 

* PROGRAM SACM1.C vsn 1.0 

* WRITTEN BY: Steven Crocker 

* LAST UPDATE: August 6, 1991 

x 

* This program takes SACM data from the Heard Island West Coast 

* Array (upper instrument package) and condenses it. 

x 

* Pressure 1s converted to depth. 

x 

x =The conductivity is not processed. All values output are 60 

* second averages of the input. 

x 

tbe / 

#define VERSION ANSI /* either ANSI or UNIX x/ 
Hdefine SELECT OUTPUT OFF /* ON or OFF x / 
Hdefine C 25 /k max length for file names */ 
H#define imBUFFER 1200 /* input buffer size x/ 
H#define outBUFFER 50 /* output buffer size x/ 


H#define PI 3.14159265359 
Hdefine RADIAN 57.2957795131 


Hinclude<stdio.h> 
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Hinclude <malloc.h> 
#include<math.h> 


#if (VERSION==ANS1) 
#include<stdlib.h> 

void ExitOnError(char error_txtl]); 
#endif 


main() 


t 


FILE *sacm_fp, *out_fp; 
char c, sacm_file{C], out_filelC]; 
int i, j, count, loop_count, inBUFFER full, out_flag([6]; 


float t1, avgTilt, avgPress, avgVN, avgVE, avgTemp, avgAngle, 
avgVelocity, junk; 


static float inTiltCinBUFFER], inPressCinBUFFER], inVNCinBUFFERI, 
INVECInBUFFER], inTimeCLinBUFFER], inTempCinBUFFERI, 
outTimeCoutBUFFER], outTiltCoutBUFFER], 
outPressCoutBUFFER], outAnglefoutBUFFER], 
outVelocityLoutBUFFER], outTempCoutBUFFER]; 


printf("\n\n\nEnter file name for SACM data. \n"); 
scanf("Zs",sacm_file); 

printf ce’ \n\n'; 

printf("Enter file name for output.\n"); 
scanf("Z%s",out_file); 

printf (C"\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n"); 
printf("The following pertains to output file format.\n\n"); 


printf("Units of measure: Time index is minutes. \n"); 
printt Angles in degrees\n"); 
printf(" Compass direction (TRUE)\n"); 
printf(" Velocity in meters/second\n"); 
printf(" Temperature in degrees C\n"); 
Brine? <~ Depth in meters\n\n'"); 


printf("NOTE: Output column order is same as above.\n\n\n"); 
for (1=0;1<6;1++) out_flagliJ=1; 


#if (SELECT _OUTPUT==0N) 
primtt("inelude TIME? (y orn \n..- 
scanf("Z%s",8&c); 

if (c==89 || c==121) out_flag([0]=1; 
else out_flagl0J=0; 

printf("\n\n"); 


print? © include TILT? <y or n)\no, 
scanf("Z%s",8&c); 

if(c==89 |; c==121) out_flagl1]=1; 
else out_flag{1J=0; 

printf ("\n\n"); 


printf("Include CURRENT DIRECTION? (y or n)\n"); 
scanf("Z%s",&c); 

if(c==89 |; c==121) out_flagl2]=1; 

else out_flag{2]=0; 

print? <*\n\ni-); 


printf("Include CURRENT VELOCITY? (y or n)\n"); 
scanf("%s",&c); 

if(c==89 ;; c==121) out_flagl3]=1; 

else out_flagf{3]=0; 

print? \nNn.)- 


printf("Include TEMPERATURE? (y or n)\n"); 
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scanf("%s" ,&c); 

if(c==89 {; c==121) out_flagl4]=1; 
else out_flagl4]=0; 

print? (\n\n"); 


printf("Include DEPTH? (Cy or n)\n"); 
scanf("%s",&c); 

if(c==89 || c==121) out_flagl(5]=1; 
else out_flagl5J=0; 

#endif 


INBUFFER_ful l=inBUFFER; 

Loop_count=0; 
if ((sacm_fp=fopen(sacm_file, "rt")) != NULL) ; 
else ExitOnError("Error opening input data file"); 
if (Cout_fp=fopen(out_file, "wt")) != NULL) ; 
else ExitOnError("Error opening output data file"); 


while (1) 
{ 
for (i=0; 1<inBUFFER; i++) 
{ 
if (fscanf(sacm_fp, "“f “Ff “fF “fF ~F ZF AFAN", 
SinTimelil,&inVNCiJ, SinVECiI, &inTemplid, 
SintiltCil, &inPressCil, &junk) != EOF) ; 
else inBUFFER_full=i-2; 
} /* for */ 
}=0; 
for (1=0;1<iMBUFFER_full; i++) 
< 
t1=inTimelid; 
avgTilt=0.0; 
avgPress=0.0; 
avgVN=0.0; 
avgVE=0.0; 
avgTemp=0.0; 
count=0; 
while(CinTimeLi]<t1+60.0) && (1<iIMBUFFER_full)) 
t 
avgTilt=avgTilt + inTiltlij/10.0; 
avgPress=avgPress + inPressliJ*1000.0; 
avgVN=avgVN + inVNC1)/1000.0; 
avgVE=avagVE + inVECLi]/1000.0; 
avgTemp=avgTemp + inTempCij]/100.0; 
J++; 
count+t+; 
} /* while */ 
avagTilt=avgTilt/( float) count; 
avgPress=avgPress/(float) count ; 
avgTemp=avgTemp/ (float) count; 
avgVN=avaVN/ (float) count; 
avgVE=avgVE/( float) count; 
avgVelocity=sqrt (avgVN*avgVNtavgVE*avaVE) ; 
avgAngle=atane(avgVE, avgVN)*XRADIAN; 
if CavgAngle<=0.0) avgAngle=360.0+avgAngle; 
out Timel jJ=(float) (j+1)+outBUFFER* loop count; 
outTiltljJ=avgTilt; 
outPress[j ]=(avgPress-101325.0)/(9.80665*1026.0); 
out TempL j J=avgT emp; 
outVelocityL[jl=avgVelocity; 
outAngle[jJ=avgdngle; 
ian 
} /* for */ 
if CinBUFFER full != inBUTFER) j--; 
for (1=0;1<j; i++) 
< 
if Cout_flag(O]) fprintf(out_fp, "%6.0f", outTimeli]); 


TD 


} 


if Cout_flagl1]) fprintt(out fp, 7 esc. otmeout litte i), 
if (out_flagl2]) fprintf(out_fp, " %8.3f", outAngleli]); 


if (out_flagf[3]) fprintf(out_fp, " %8.6f", outVelocitylil); 
if (Cout_flag(4]) fprintfCout_fp, "  %8.5f", outTemplil); 
if (out_flag[5]) fprintfCout_fp, " %8.0f", outPressliJ); 
Tprinttcout tp, vn, 

Fs k Or */ 


if CinmBUFFER full != inBUFFER) exit(0); 
loop_count+t+; 

} /* while */ 

exitc))- 


Hif (VERSION==ANSI) 

void ExitOnError(char error_txtl]) 
Helif (VERSION==UNIX) 
ExitOnError(error_txt) 

char error_txtCI; 

Hendif 


{ 


De 


fprintf(stderr,"Program run-time error ...\n"); 
fprintf(stderr,"%s\n",error_txt); 
fprintf(stderr,"...now exiting to system...\n"); 
exit(1); 


SACM2.C 


This program takes the output of sacml.c as input. ie 


locates and copies a user defined subset of the tilt data for 


use in beamforming. 


[kkk 

x PROGRAM SACM2.C vsn 1.0 

x WRITTEN BY: Steven Crocker 

* LAST UPDATE: August 6, 1991 

x 

x This program takes data processed by SACM1.C and cuts 

x a user defined segment from it. The segment is retained 

x as the output data file. The input file is not affected. 

x The output time base may either be normalized to 1 or may 

x retain the original values. 

x 

x The option is given to accept a default output format. This 

x output format coincides with the required input format to the 
x time domain modal beamformer. 

x 

x The output elements are selectable if desired. 

KKKKK / 

Hdefine VERSION ANSI /k either ANSI or UNIX x/ 
H#define C 25 /x max length for file names */ 


#include<stdio.h> 
#Hinclude <malloc.h> 


Hif (VERSION==ANSI) 
#include<stdlib.h> 

void ExitOnError(char error_txtClJ); 
Hendif 
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main() 
{ 
FILE xinput_fp, *out_fp; 
char c, Default, input_file[C], out_filelC]; 
int 1=0, out1, outN, out_flag[7J; 
float Tilt, Depth, Time, Temp, Angle, Velocity; 


printf("\n\n\nEnter file name for input data.\n"); 
scanf("%s", input_file); 

if (Cinput_fp=fopen(input_file, "rt")) != NULL) ; 

else ExitOnError("Error opening input file"); 

Brimtt<(’ \n\n"); 


printf("Enter file name for output. \n"); 
scanf("%s",out_file); 
if (Cout_fp=fopen(out_file, "wt")) != NULL) ; 
else ExitOnError("Error opening output file"); 


printf("\n\n\nDo you want the default output configuration?\n"); 
printf("This option formats the output for use in the beamformer. \n"); 
scanf("%s" ,Default); 
printf¢"\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n") ; 


if(Default '= ‘n' && Default != 'N') 

{ 
printf("The following pertains to output file format.\n\n"); 
printf("Units of measure: Time Index is minutes. \n"); 
primer C.. Angles in degrees\n"); 
print th | Compass direction (TRUE)\n"); 
printtC’ Velocity in meters/second\n") ; 
primttc: Temperature in degrees C\n"); 
printf(" Depth in meters\n\n"); 


printf("NOTE: Output column order is same as above. \n\n\n"); 


printf("Include TIME? Cy or n)\n"); 
scanf("%s",¢); 

if(c==89 |; c==121) out_flagl[0J=1; 
else out_flag[0J=0; 

Ph imiet (© \mAn: ):; 


printf("Include TILT? (¢y or n)\n"); 
scant G'As  c): 

if(c==89 |; c==121) out_flagl1]=1; 
else out_flag[1]=0; 

Peiigtic \n\ni-); 


printf("Include CURRENT DIRECTION? (y or n)\n"); 
Scant« AS.,¢c)> 

if (c==89 {| c==121) out_flagl[2]=1; 

else out_flag[2]=0; 

Pram ne \nAm:)); 


printf("Include CURRENT SPEED? (y or n)\n"); 
scant("4s",¢); 

if(c==89 |; c==121) out_flagl[3]=1; 

else out_flag[3]=0; 

peintr ¢’ \n\n'); 


printf("Include TEMPERATURE? (y or n)\n"); 
scanf("%s",c); 

if(c==89 |; c==121) out_flag[4]=1; 

else out_flag[4]=0; 

Printi< \a\n')- 


printf("Include DEPTH? (y or n)\n"); 
scanf("%s",c); 
if(c==89 |; c==121) out_flag[5]=1; 
else out_flag[5]=0; 

} /* if */ 


a 


else 


out_flagl0J=0; 

out_flagl1J=1; 

out_flagl2)=1; 

out_flagl3J=0; 

out_flagl4]=0; 

out_flag(5J=1; 
+ /* else */ 


printf("\n\nWhat time index of the input file should be the first\n"); 
printf ("element of the output file?\n"); 

scanf("%i",&out1); 

Brintt © \Oin, 


printf("How many elements should the output file contain?\n"); 
scanf ("%i",&outN) ; 
printf("\n\n"); 


if (Default != 'n' && Default != 'N') 


{ 
printf ("Should the time base be normalized to begin at t=1?\n"); 
scanf("%s",c); 
if(c==89 |; c==121) out_flag(6J=1; 
else out_flagl6J=0; 
Print? (\n\n\n 
FS * 11a 
while (1) 
. 


if (fscanfCinput_ fp, "%f ~2f “fF “AF ZF ZF \n", 
&Time, &Tilt, &Angle, &Velocity, &Temp, &Depth) != EOF) 


{ 
if((Time >= (float)outl) && (Time < (float) Cout1+outN))) 
{ 
ifCout_flag(O] && !out_flagl6]) 
fprintt (OULSTp, -41. >) lime), 
else if Cout_flaglO] && out_flagl[6]) 
fprintt(out_ tp, 41- (fT loaw Citi), 
if Cout_tlagl1)) forintf outlier 
if Cout_flagl2]) fprintt(out_tp,) 4 <1 7 Angle), 
if Cout_flagl3J) fprintf(out_fp, " <%f", Velocity); 
if <out_flagl(4J) fprintf(out_fp, " %f", Temp); 
it Cout flaglis)) torinttcoutstp, “ es 41vepEn, 
fprinttCout fp) ene), 
i++ > 
Feet ey 
else if(Time >= (float) Cout1+outN)) 
{ 
fclose(input_fp); 
fclose(out_fp); 
exit(Q); 
+ /* else if */ 
>} /k if */ 
else 
a 


fclose(input_fp); 
fclose(out_fp); 
ExitOnError("EOF encountered in input data file"); 
} /*® else */ 
} /k while */ 
exit(1); 
} 


#if (VERSION==ANS]) 

void ExitOnError(char error_txt(J) 
Helif (VERSION==UNIX) 
ExitOnError(error_txt) 
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char error_txt(]; 


Hendif 

{ 
fprintf(stderr,"Program run-time error ...\n"); 
fprintf(stderr, "%s\n",error_txt); 
fprintf(stderr,"...now exiting to system...\n"); 
exit(1); 

} 


E. ARRAYTEST.C 


This program performs various statistical tests on each 


element of the acoustic array. It is useful for isolating 


hydrophones which have failed during the experiment. The 


statistical tests are an adaptation of those 


Reference ll. 


| xkekx 

* PROGRAM: ARRAYTST vsn 1.0 

* WRITTEN BY: Steven Crocker 

* LAST UPDATE: October 21, 1991 


x 

* 

Keke / 

Hdefine UNIX VERSION /*® either ANSI or UNIX 

Hdefine CHANNELS 32 /* number of hydrophones on array 
Hdefine F_SAMPLE 228 /* sampling frequency 

#define BUFFER TIME 60 /* duration of time averages 


Hinclude<stdio.h> 
Hinclude<malloc.h> 
#include<math.h> 


Hif defined ( ANSI ) 
Hinclude<float.h> 
#Hinclude<stdlib.h> 
int getInput (void); 
int moments(float *data, int n, float *ave, float *ave2, float *adev, 
float *sdev, float *svar, float *skew, float *curt) 

float *vector(int length); 
float **matrix(int row, int col); 
int free_matrix(float **m, int row); 
void ExitOnError(char error_txtl]); 
Helif defined ( UNIX ) 
int vector(), 

matrix(), 

getInput(), 

moments(), 

free matrix, 

ExitOnError(); 
Hendif 


/* Global Variable Declarations */ 
float **inSOUND; 
int firstBuffer=1; 


maint) 


{ 
int i, j=1, k, flagl8); 
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x/ 
*/ 
*/ 
*/ 


found in 


float ave, ave2, adev, sdev, svar, skew, curt, *temp; 
char fileName{80], key; 
FILE *fp1, *fp2, *fp3, *fp4, *fp5, *fp6, *fp7; 


inSOUND=(float **)matrix(CHANNELS, BUFFER _TIME*F_SAMPLE); 
temp=(float *)vector(BUFFER_TIMEXF_SAMPLE) ; 


for (1=0;1<=7; i++) flagliJ=0; 
printf("Answer y or n to the following output options. \n"); 
printf ("Average value ( over %i seconds ): ", BUFFER_TIME); 
scanf("%s", &key); 
if (key==89 {| key==121) 
{ 
flagl1IJ=1; 
printf("Enter file name: "); 
scanf("%s", fileName); 
PrINnGT Gavan, 2; 
if ((fpl=fopen(fileName, "wt")) == NULL) 
ExitOnError("Error opening file."); 
} /* if */ 


printf("Mean squared value ( over %i seconds ): ", BUFFER_TIME); 
scanf("Z%s", &key); 
if (key==89 || key==121) 
< 
flagl23=1; 
printf("Enter file name: "); 
scanf("Z%s", fileName) ; 
DEINtiC. \A\M, 2; 
if ((fp2=fopen(fileName, "wt'")) == NULL) 
ExitOnError("Error opening file."); 
>} /* if */ 


printf("Average deviation ( over %i seconds ): ", BUFFER_TIME); 
scanf("Z%s", &key); 
if (key==89 |; key==121) 
{ 
f lag{3J=1; 
printf("Enter file name: "); 
scanf("%s", fileName) ; 
Dr IM ta | nN 
if ((fp3=fopen(fileName, "“wt")) == NULL) 
ExitOnError("Error opening file."); 
} /k if */ 


printf("Standard deviation ( over %i seconds ): ", BUFFER_TIME); 
scanf("%s", &key); 
if (key==89 |; key==121) 
< 
flag(4J=1; 
printf("Enter file name: "); 
scanf("Z%s", fileName); 
Orin cL > NM Nee 
if ((fp4=fopen(fileName, "wt")) == NULL) 
ExitOnError("Error opening file."); 
yo Roa ies 


printf("Variance ( over %i seconds ): ", BUFFER TIME); 
scanf("Z%s", &key); 
if (key==89 |; key==121) 
{ 
flagl53=1; 
printf("Enter file name: "); 
scanf("%s", fileName); 
PrInefe: \n\n')- 
if ((Cfp5=fopen(fileName, "wt")) == NULL) 
ExitOnError("Error opening file."); 
} /k if */ 
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printf("Skewness ( over %i seconds ): ", BUFFER_TIME); 
scanf("%s"", &key); 
if (key==89 {| key==121) 
{ 
flagl6)=1; 
printf("Enter file name: "); 
scanf("%s", fileName); 
printf C"\n\n"); 
if ((fp6=fopen(fileName, "wt")) == NULL) 
ExitOnError("Error opening file."); 
} /k if */ 


printf("Kurtosis ( over %i seconds ): ", BUFFER_TIME); 
scanf("%s", &key); 
if (key==89 |; key==121) 
{ 
flag(7J=1; 
printf("Enter file name: "); 
scanf("%s", fileName); 
ramen c \n\n); 
if ((fp7=fopen(fileName, "wt")) == NULL) 
ExitOnError("Error opening file."); 
} /* if */ 


while (}<120) 

{ 
get Input (); 
mectlagtti) fprintt(ip1,"%1", 3); 
MmarGtlageed)) fprintt({pe,"41",3); 
if (flag(3]) fprintf(fp3,"%i",j); 
fect lagiad)) tprintt (tps, "41", 3); 
mictlag(51) fprintt(fp5,"%1",j); 
if (flagl6]) fprintf(fp6,"%i",j); 
miect lagi7]) fprintf(fp7,"%41i",)); 


for (1=1; 1<=CHANNELS; i++) 
c 
for(k=1;k<=BUFFER_TIMEXF_SAMPLE; k++) 
tempCkJ=inSOUNDCiICkI; 


moments (temp,BUFFER_TIMEXF_SAMPLE, &ave,&avec, 
Radev,&sdev,&svar,&skew,&curt); 

am (flaglti) forintf(fpi1, " “<f ", ave); 

if (flagl2]) fprintf(fp2, " “~f ", ave2); 

iMectlacksa atprintt(tps, “ 4f ", adev); 

if (flagl4]) fprintf(fp4, " %f ", sdev); 

if ct laglod) tprinti(tps, ~ 4T ", svar); 

if (flagl(6]) fprintf(fpé6, " %f ", skew); 

iiectlagerap fordatt ipr. © 7f “,curt); 

+} /* for */ 


itieccl doe le tprintt (ipipJ\n).: 
if (flagl2]) fprintf(fp2,"\n"); 
imeect laglh3))) forintt (fp3 ,7\nl; 
if (flagl4J) fprintf(fp4,"\n"); 
it vaGk dU) aprintt (fps, "\n"); 
if (flagl6]) fprintf(fp6,"\n"); 
imecteagi/?)) fprintt(to7,"\n"); 
firstBuffer=0; 

jtt; 

} /* while */ 


[oC ddd END main dodging ook / 
[IKARIA IKKE IKAIAAAA AAA IAAI AAAI AIA AKIRA AAI AAAI IAAAN / 
[xkxkk 
* FUNCTION: getInput() 
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This function handles all acoustic input. It also provides one of 
two normal process terminations available in the program. 


x 

x 

x 

x 

x Arguments: none 

x 

x Return value: 0 

x 

x Functions called: vector() ExitOnError() 
x 

x Definitions called: ANSI UNIX 

. F_SAMPLE CHANNELS 

* BUF FER_TIME 

x 

* Global variables called: insOUNDCILI firstBuffer 
x 

x Significant memory allocation:  diskBuffer{] 

*x 

tks / 


Hif defined (¢ ANSI ) 
int getInput (void) 
Helif defined ( UNIX ) 
get Input() 
#endif 
{ 
int i, j, items, buffer; 
float *diskBuffer; 
char fileName[80); 
Static FILE *fpStatic; 


if (firstBuffer) 
{ 
printf("Enter file name for input acoustic data: "); 


if((scanf("%s", fileName) )==EOF) 
ExitOnError("Fatal error in scanf()"); 


Drancorc NAAN 2 


if (CfpStatic=fopen(fileName, "rb")) == NULL) 
ExitOnError("Error opening INPUT ACOUSTIC data file."); 
> f* ates 


/* Memory Allocation */ 
buffer=BUFFER_TIMEXF_SAMPLE*CHANNELS; 
diskBuf fer=(float*)vector(buf fer); 


items=fread((char*) (diskBuffer+1), sizeof(float), buffer, fpStatic); 
if Citems==buffer) /*continue*/; 
else if(ferror(fpStatic) != 0) 
ExitOnError("Error encountered while reading input acoustic data"); 
else if(feof(fpStatic) != 0) 
< 
PP INTE CUA CRRRKREKAAEKK RAK KAKKIEIKRKIAAKK IDI IISIAIAIAAAAE IIR 1) 5 
Print tc. \e End of File reached: EXECUTION COMPLETE\n"); 
DP TNCE COCKER IRI R IAAI IA KI IAAI III III DADA IAI IAI IIIA AAAI IAI AIAN PH) 
fclose(fpStatic); 
exit (0); 
+} /* else if */ 
else ExitOnError("Unknown error handling acoustic input file."); 


for (1=1; 1<=BUFFER_TIMEXF_SAMPLE; i++) 
{ 
for (j=1; j<=CHANNELS; j++) 
inSOUND(LjICiJ = disk3ufferlCHANNELS*(1-1)+)jJ; 
> 1% for %/ 


/* Deallocate Memory */ 
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free((char*)diskBuf fer); 
return( 0 ); 


> iti itt itis / 


/kkkkk END getInput **xkk/ 
[aid / 


/Kxkkr 

* FUNCTION: moments() 

* 

* This function calculates the mean, mean squared value, average 
* deviation, standard deviation, variance, skewness and kurtosis 
* of a data vector. 

* 

* Arguments: data n 

* ave ave2 
* adev sdev 
x svar skew 
* curt 

& 

* Return value: 0 

* 

* Functions called: ExitOnError() 

* 

* Definitions called: ANSI ; UNIX 
* 

* Global variables called: none 

* 

x Significant memory allocation: none 

* 

tot / 


#if defined ( ANSI ) 
int moments(float *data, int n, float *ave, float *ave2, float *adev, 


float xsdev, float *svar, float *skew, float *curt) 


Helif defined ( UNIX ) 
moments(data,n,ave,ave2,adev, sdev, svar,skew, curt) 


mt n; 
float *data, *ave, x*ave2, *adev, xsdev, *svar, xskew, *curt; 
#Hendif 
{ 
Met); 
float s,p; 
if (n<=1) ExitOnError("n must be at least 2 in moment()"); 
s=0.0; 
kave2=0.0; 


for (j=1; j)<=n; jt+) 
{ 


s += datal)J; 

kave2d += dataljl*datalj]; 
} /* for */ 
kave=s/n; 
kave2 /= n; 
kadev=(xsvar)=(*skew)=(xcurt)=0.0; 
for (j=1; j<=n; j++) 
{ 


kadev += fabs(s=data[jJ-(*ave)); 
ksvar += (p=s%*s); 
xskew += (p *= s); 
kcurt += (p *= s); 
} /* for */ 
kadevy /= n; 
ksvar /= (n-1); 
ksdev=sqrt (*svar); 
if (*svar) 
q 
kskew /= (nk(xsvar)*(x%sdev)); 
kcurt=(Xcurt) /(nk(xsvar)* (xsvar) )-3.0; 
> /* if */ 
else ExitOnError("No skew/kurtosis when variance = O (in moment())"); 


oe 


return( 0 ); 
3 Dodo oI / 


[/kkkkk END moment ***k*/ 
LEKKKKKEKKKEEEKEEEEEEES / 


| kikkk 

* FUNCTION: vector() 

; This function allocates memory for UNIT OFFSET vectors. 
: Arguments: Length 

Return value: ky 

* Functions called: ExitOnError() 
Definitions called: ANSI UNIX 
* Global variables called: none 

Significant memory allocation: vCJ 


tke / 
Hif defined ( ANSI ) 
float *vector(int length) 
Helif defined ( UNIX ) 
vector (Length) 
int length; 
#Hendif 
< 

float xv; 


if ((v=(float*)malloc((length+1)*sizeof (float) ))==NULL) 
ExitOnError("Memory allocation failure in vector()."); 
#if defined ( ANSI ) 
return v; 
Helif defined ( UNIX ) 
return (long int)v; 
Hendi f 
FY [KARRI TID IAIKAEHIAKK / 
[/xkkkkK END vector *xkkx/ 
[KKK KKKIA KIARA IA II KK / 
[kkKKK 


: FUNCTION: matrix() 

; This function allocates memory for UNIT OFFSET 2-D arrays. 
A Arguments: row col 

* Return value: *km 

: Functions called: ExitOnError() 

Definitions called: ANSI UNIX 
Global variables called: none 

x Significant memory allocation: m{IJCJ 

KeeKK/ 


#if defined ( ANSI ) 
float *xmatrix(int row, int col) 
Helif defined ( UNIX ) 
matrix (row, col) 
int row, col; 
Hendif 
{ 
Nt a 
float **m; 
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if ((m=(float**)mal loc ( (unsigned) (row+1)*sizeof (float*)))==NULL) 
ExitOnError("Allocation failure 1 in matrix()"); 
for (1=1; 1<=row; 1++) 
a 
if ((mLiJ=(float*)mal loc( (unsigned) (col+1)*sizeof (float) ) )==NULL) 
ExitOnError("Allocation failure 2 in matrix()"); 
} /* for */ 
Hif defined ( ANSI ) 
return m; 
Helif defined ( UNIX ) 
return (long int)m; 
Hendif 
3 Dito / 
/kkkkk END matrix *eekk/ 
[oid / 


[xkxkk 

id FUNCTION: free_matrix() 

* This function deallocates memory from UNIT OFFSET 2-D arrays. 
; Arguments: mC ICI row 
Return value: 0 

* Functions called: none 

Definitions called: ANSI UNIX 

: Global variables called: none 

* Significant memory allocation: none 

cht 


Hif defined ( ANSI ) 

void free_matrix(float **m, int row) 
Helif defined ( UNIX ) 
free_matrix(m, row) 

float **m; 

Int row; 

H#endif 

. 


int 4; 


for(i=row; 1>=1; 1--) 
free((char*)m(1]); 
free((charx)m); 
return( 0 ); 
9) I IIIT IATA IAI / 
/kkkkk END free matrix *kkkk/ 
[SKAAKIKKEKEEKIK ERIK AKI / 


| kkkkk 

A FUNCTION: ExitOnError() 

; This function performs an abnormal process termination. 
F Arguments: error _txtC] 

* =Return value: none 

* 

‘ Functions called: none 

* Definitions called: ANSI UNIX 
: Global variables called: none 

* 

* Significant memory allocation: none 

reek / 


85 


#if defined ( ANSI ) 

void ExitOnError(char error_txt{]) 
Helif defined ( UNIX ) 
ExitOnError(error_txt) 

char error _txt(); 


fendif 

c 
fprintf(stderr,"Program run-time error ...\n"); 
fprintf(stderr, "%s\n",error_txt); 
fprintf(stderr,"...now exiting to system...\n"); 
exit (0); 


} ooo obi / 


/kekkKK END ExitOnError *xkkk/ 
[on oo oo / 
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